Rok 2018 zbliża się powoli ku końcowi. Jest to często czas podsumowań i refleksji. Pomyślałem więc, że również pozwolę sobie na małą refleksję odnośnie do jakości powietrza w Polsce w roku 2017 (sic!).
Dlaczego rok 2017? Bo akurat w momencie pisana tego artykułu jest to najnowsza paczka danych zawierająca archiwalne pomiary dostępna na stronie GIOŚ. Fakt, możemy dobrać się do aktualnych danych GIOŚ przez API, ale pozwala ono na dostanie się do ostatnich trzech dni kalendarzowych. Jeśli więc nie zbieraliśmy danych wcześniej w 2018 roku, to żeby zrobić podsumowanie roku 2018, będziemy czekać na odpowiednią paczkę danych, która pojawi się jakoś w trakcie roku 2019. Proste.
Ideą podsumowania roku 2017 będzie zbudowanie listy miast i ułożenie jej według jakości powietrza. Coś podobnego robiłem już wcześniej, tym razem będą to nieco inne dane, więcej miast i trochę dodatkowych przemyśleń. Do dzieła więc.
Dane
Dane, które będę wykorzystywał, pochodzą ze wspomnianego powyżej archiwum. Znajdziemy je w pliku 2017.zip. Potrzebny nam też będzie plik z informacją o współrzędnych stacji pomiarowych – Metadane_wer20180829.xlsx. Do pełni szczęścia potrzebny nam jeszcze plik z listą miast i ich współrzędnymi. Plik ten przygotowałem przy pomocy Wikipedii oraz OpenStreetMaps, a pobierzesz go tutaj (po pobraniu zmień rozszerzenie z txt na csv).
Proces
Proces tworzenia podsumowania składa się z kilku etapów. Pierwszy etap to wczytanie i przetworzenie danych. Tutaj wszystko oparłem o pętlę, która przechodzi po każdym z czynników szkodzących. W takiej pętli mam kolejną pętlę, która przechodzi po liście miast. Dla każdego miasta obliczam odległość do wszystkich stacji mierzących ten czynnik szkodzący i wybieram najbliższą. Wtedy pobieram dane dla tego czynnika z tej stacji i wyliczam wartość średnią z całego roku. Tak uzyskaną średnią zapisuję w ramce danych.
Jak można się domyślić 7 czynników, 930 miast i około 100 stacji na każdy czynnik daje trochę złożoności obliczeniowej. I to jest na razie, chyba najsłabszy element mojej analizy. Mam już pewien pomysł, jak przyspieszyć ten element, ale o tym na końcu.
Drugim elementem jest normalizacja odczytów. Jako że czynniki szkodzące mierzone są w różnych jednostkach i różne ich stężenie ma odmienny wpływ na zdrowie człowieka, postanowiłem je znormalizować, żeby traktować je porównywalnie. Wartość „bardzo zły” da po normalizacji wartość 1.
Ostatnim elementem jest wyznaczenie wartości średniej tak znormalizowanych czynników i posortowanie tabeli. Proste, nie?
Problem
Otóż nie jest to takie proste. Fakt, aktualnie zadowalam się powyższym rozwiązaniem. Ale to dlatego, że nie wymyśliłem jeszcze jak rozwiązać pewien problem. Problemem tym jest odległość miasta od czujnika, a raczej problem braku uwzględnienia tej odległości.
Ale o co chodzi? Otóż okazuje się, że GIOŚ dostarcza całkiem mało niektórych rodzajów czujników. Najbardziej „lubiany” w Polsce czynnik PM2.5 jest monitorowany w około 50 miejscach w Polsce. Nie udało mi się ustalić dlaczego. Skoro mamy 930 miast w Polsce, a czujników 50, to istnieje całkiem spora szansa, że jakieś miasta będą położone dość daleko od tego typu czujników. Myślę, że warto byłoby uwzględnić ten „efekt” w podsumowaniu, ale jak?
Pomysłem, który mam, jest wyznaczenie średniej ważonej, gdzie wagi są odległością. Im mniejsza odległość, tym waga jest większa. Jako hipoteza brzmi to dobrze i nawet pokusiłem się o sprawdzenie rankingu z uwzględnieniem wag znormalizowanych po odległości wszystkich czujników i alternatywnie znormalizowanych po odległościach pogrupowanych czujników. Oczywiście Python przyjmie wszystko, więc taki ranking powstał. Co się okazało, jako jedne z najlepszych powychodziły miasta, które miały bardzo daleko do czujników PM2.5. Dlaczego? No bo PM2.5 jest w Polsce najbardziej problematyczne. Pozostałe czynniki nie są aż takie złe. No i mamy więcej czujników, które je monitorują. Do średniej ważonej wpadły więc niskie wartości (bo np. nie mamy problemu z O3) z dużymi wagami (bo są blisko), a PM2.5 mimo iż duże, miało małą wagę, bo jest mierzone daleko. Gdzie tkwi więc haczyk?
Jeśli handlowałbym np. nieruchomościami i byłbym na tyle dużym graczem, że sterowałbym rynkiem, to moja metryka (średnia ważona) mogłaby wytworzyć popyt na mieszkania. W interesie takich miast byłoby utrzymanie wysokiej pozycji w rankingu najlepszego powietrza, bo wskazywałbym je jako dobre miejsce do osiedlenia się. Nie inwestowałyby więc w czujniki PM2.5, bo skoro, w rankingu jest dobrze, to jeśli okaże się, że w mieście jest źle, to waga tego pomiaru by wzrosła. Wtedy spadłyby w rankingu i by straciły. Taka metryka ze średnią ważoną jest więc w mojej opinii w tej chwili zła, bo optymalizacja zakłada nieinwestowanie w czujniki najbardziej problematycznego czynnika.
Nie wymyśliłem jeszcze jak rozgryźć ten problem. Nie chciałbym ustawiać żadnych sztywnych granic, progów i tym podobnych. W idealnym świecie każde miasto powinno mieć w geometrycznym centrum baterię 7 czujników. Problem odległości przestałby wtedy istnieć. I tak pewnie będzie za jakiś czas. W tej chwili myślę, że warto byłoby jakoś uwzględnić tę odległość. Może Ty masz jakiś pomysł jak to ugryźć?
Wyniki
Zacznijmy może od odległości od czujników. Gdy wyszukiwałem najbliższy czujnik danego czynnika dla każdego miasta, oprócz zapisania sobie jego średniego odczytu z całego roku, zapisywałem również jego odległość. Mogę wiec sprawdzić, jak ma się np. 5 największych odległości dla każdego czynnika:
PM25 PM25_dist
city
Mielno (ZP) 28.3591 126.39
Koszalin (ZP) 28.3591 117.773
Tychowo (ZP) 28.3591 114.676
Barwice (ZP) 28.3591 114.397
Białogard (ZP) 20.0019 111.765
PM10 PM10_dist
city
Terespol (LB) 29.0385 91.8119
Przasnysz (MZ) 24.9663 80.6835
Włodawa (LB) 32.5306 75.6202
Ciechanów (MZ) 31.5396 71.6462
Mława (MZ) 26.855 70.4742
C6H6 C6H6_dist
city
Krynki (PL) 1.52227 220.139
Sokółka (PL) 1.03138 204.031
Michałowo (PL) 1.52227 198.077
Czarna Białostocka (PL) 1.03138 193.413
Lipsk (PL) 1.03138 192.507
NO2 NO2_dist
city
Sokołów Podlaski (MZ) 11.595 74.4007
Łobez (ZP) 14.5201 69.7483
Gryfice (ZP) 22.9793 68.5706
Dziwnów (ZP) 22.9793 67.2836
Resko (ZP) 14.5201 65.7629
SO2 SO2_dist
city
Sokołów Podlaski (MZ) 5.3082 74.4007
Łobez (ZP) 2.3102 69.7483
Gryfice (ZP) 4.67442 68.5706
Dziwnów (ZP) 4.67442 67.2836
Resko (ZP) 2.3102 65.7629
O3 O3_dist
city
Czaplinek (ZP) 44.5335 106.062
Kołobrzeg (ZP) 50.2521 100.889
Połczyn-Zdrój (ZP) 50.2521 99.5969
Gościno (ZP) 44.5335 99.2115
Wałcz (ZP) 42.6561 98.4846
CO CO_dist
city
Terespol (LB) 0.403192 116.04
Łaszczów (LB) 0.244055 114.43
Tyszowce (LB) 0.403192 106.784
Hrubieszów (LB) 0.403192 105.403
Ostrołęka (MZ) 0.432785 102.196
W powyższej tabeli widzimy, że odległości 60+ kilometrów nie są niczym zaskakujące. Najbardziej zastanawiają liczby powyżej 100, a nawet 200 kilometrów. Okazuje się zatem (o ile nie popełniłem żadnego błędu), że mieszkańcy miasta Krynki w 2017 roku najbliższy państwowy czujnik benzenu mieli w odległości około 220 kilometrów. Z tego, co widziałem w aktualnych danych, to Białystok ma już taki czujnik, więc sytuacja wygląda nieco lepiej. Ale o tym w kolejnym artykule – poszukiwanie czujnikowych białych plam.
Okej, przejdźmy więc do najczystszych i najbrudniejszych miast w Polsce. Przypominam metodologię – każde ma przyporządkowaną grupę najbliższych czujników – po jednym na każdy czynnik szkodliwy. Bierzemy średnią z cogodzinnych średnich z całego roku 2017. Wartości te normalizujemy tak, że granica „bardzo zły” daje jeden. Na końcu wyliczamy średnią tych znormalizowanych wartości i sortujemy tabelę.
Najczystsze miasta w Polsce w 2017:
Najbrudniejsze miasta w Polsce w 2017:
Jak widać, województwo pomorskie całkiem dobrze przędzie w jakości powietrza. Cóż, wynika to pewnie z faktu, że województwo to przylega do morza. Możemy więc sobie wyobrazić, że nawet jeśli mieszkańcy tych miast nieźle ładowaliby do kotłów, to i tak morze stanowi całkiem spory rezerwuar pustej przestrzeni, która może przyjąć te czynnik szkodliwe. Nie mówiąc już o tym, że samo morze często powoduje dość sporą intensywność wiatrów. Może o to więc chodzi?
Co do najbrudniejszych miast cóż, nie jestem specjalistą w tej dziedzinie, ale pamiętam, że uczono mnie w szkole, że właśnie w tamtych obszarach leżących na terenie województw małopolskiego i śląskiego mamy w Polsce najwięcej przedsiębiorstw wydobywczych, przetwórczych i produkcyjnych. Oczywiście chodzi mi o przemysł ciężki i surowce geologiczne. Czyli po prostu kopalnie, huty i fabryki. Plus, z jednej strony obszar ten jest „ograniczony” przez Karpaty, więc to też może niejako „blokować” rozchodzenie się smogu. Ale to tylko takie zgadywanki, nie traktuj więc tego i powyższego akapitu na poważnie ;).
Cały proces wyznaczania miasta, jak i funkcja „debug” pokazująca informacje o czujnikach, które zostały przypisane do miasta, znajduje się w notebooku, do którego link znajdziesz na końcu artykułu. Zachęcam do samodzielnego poeksperymentowania ;).
Możemy jeszcze się przyjrzeć, z których stacji pobrane zostały dane dla poszczególnych miast. Poniższe mapy są zrzutami ekranu z notebooka z analizą.
Widzimy, że w mojej analizie największym problemem jest odległość do czujników. No bo odległość czujników od miasteczka Skała jest znacząca, nie mówiąc już o tym, że czujniki te są umieszczone w różnych obszarach. Cóż, pozostawię to więc tak jak jest, do przyszłych rozważań.
Przyszłość
Tak jak wspominałem, aktualnie najwięcej czasu marnuję na namierzenie najbliższego czujnika dla miasta. Myślę, że można by spróbować zaprząc Scikit-Learn do rozwiązania tego problemu. Na pierwszy rzut oka metoda NearestNeighbors wydaje się obiecująca. Myślę, że w jakimś kolejnym artykule spróbuję ją wykorzystać do przyspieszenia działania skryptu.
Jeżeli natrafię na jakiś fajny sposób na uwzględnienie odległości w wyznaczaniu jakości powietrza, to również postaram się go zastosować. Mam jednak przeczucie, że nie będzie miało to za dużo sensu i chyba sensowniej będzie lobbować GIOŚ o postawienie nowych czujników i o posiłkowanie się danymi z Airly (znasz może jeszcze jakieś API?).
Wciąż jeszcze nie ruszyłem modelowania czasoprzestrzennego. Kto wie, może tutaj uda mi się zbudować coś ciekawego, i być może użycie takiego modelu da nieco więcej wglądu w stan zanieczyszczenia powietrza w Polsce. Szczególnie że właśnie odkryłem, że IMGW również udostępnia historyczne paczki z danymi pogodowymi (aww yiss).
Cóż, liczę na to, że z roku na rok dostępnych będzie coraz więcej czujników. Liczę też, że nauczę się nowych technik analitycznych i będę mógł lepiej przygotować kolejne raporty. A na razie, musimy zadowolić się tym, co mamy.
Pełny kod użyty w raporcie znajduje się tutaj.
Znajdowanie najbliższych sąsiadów dla 930 miast i 100 czujników nie powinno być wyzwaniem obliczeniowym. To jest złożoność O(n*m) W kodach widzę dwa miejsca, które powodują spadek wydajności. Pierwsze to linia:
nearest_station = tmp_stations.sort_values([„dist”]).iloc[0]
Sortowanie nie jest potrzebne. Interesuje nas jeden, najbliższy czujnik, więc pewnie funkcja min wystarczy. Domyślam się, że sortowanie przydawało się w debugowaniu, do narysowania kilku najbliższych czujników. Ale przez to wydajność tej operacji spada z O(n) do O(n logn)
Drugi powód grzania procesora to linia:
dist = geopy.distance.vincenty(coords_1, coords_2).km
geopy to całkiem przydatny pakiet i bardzo precyzyjnie oblicza odległość między punktami. Ale tu nie potrzebujemy milimetrowej dokładności i uwzględniania takich niuansów jak spłaszczenie Ziemii. Przy rzadkiej sieci czujników i bardzo zgrubnych pozycjach miast nawet 10% różnica długości równoleżnikowej między południem i północą Polski nie ma znaczenia. Można Polskę zrzutować na płaszczyznę i wszystkie współrzędne geograficzne zmapować na współrzędne kartezjańskie odległości względem np. Łodzi. Dla poprawienia precyzji można nawet za mapę Polski przyjąć nie kwadrat a trapez. Odległość między punktami jest wtedy zwykłą długością odcinka pomiędzy punktami, a jej obliczenie to dwa odejmowania, dwa mnożenia i pierwiastkowanie. 930 miast * 100 czujników powinno się przeliczyć w jakąś milisekundę. sklearn.neighbours w niczym tu nie pomoże, a raczej pogorszy sprawę.
A jeśli chodzi o samo mierzenie zanieczyszczeń, to założenie o pobieraniu danych z najbliższego czujnika jest baaaaardzo daleko idące. W ten sposób Sieraków położony między jeziorami przy Puszczy Noteckiej ma poziom PM25 jak w Poznaniu. Zanieczyszczenia nie rozkładają się równomiernie, ale są skoncentrowane wokół dużych miast i ośrodków przemysłowych. A nawet lokalnie pomiary potrafią się znacząco różnić w zależności od ukształtowania terenu i pobliskich emitentów. Przy tak przyjętych założeniach wyniki naniesione na mapę pokażą plamy o jednolitej wartości skoncentrowane dookoła czujników. Podobnie mało użyteczne to podejście by było dla przedstawienie średniego wynagrodzenia, gdzie dane są dostępne tylko dla niektórych miejsc w Polsce. Tu też Sieraków miałby pewnie taką samą wartość jak Poznań i też nijak by się to nie miało do rzeczywistości.
Hej Paweł!
Dzięki za rzeczowy i konkretny komentarz.
Faktycznie całkiem zapomniałem o tym sortowaniu. Tak jak mówisz, w czasie budowania „naocznie” sobie to ogarniałem dla kilku najbliższych punktów. Przy następnej iteracji zdecydowanie to przerobię.
Co do obliczania odległości to dla mnie to pewna nowość. Szukałem trochę w internecie materiałów na ten temat i trafiłem np. na coś takiego: https://www.mkompf.com/gps/distcalc.html. I uświadomiłem sobie, że łatwo tego nie ugryzę. Dowiedziałem się m. in. tego, że zwyczajnie z twierdzenia Pitagorasa wykorzystując długość i szerokość geograficzną daleko nie zajdę. Fakt, jeśli jest tak jak mówisz ze współrzędnymi dla Polski to faktycznie nie powinno to wnosić znaczącego efektu. Na przeliczeniu wszystkiego na Łódź też nie wpadłem. Dzięki!
A co do jakości tych pomiarów to jest dokładnie tak jak mówisz, spora lipa. Dlatego też między innymi umieściłem w artykule mapki, pokazujące jak to w praktyce wygląda, bo z ramki danych średnio da się to wyczytać. Może lud się zepnie i przyciśnie trochę GIOŚ i organy nadrzędne, żeby więcej tych czujników ogarnęli? Może czujniki będą tańsze i to też pomoże. W każdym razie, wiele miast nie mierzy praktycznie niczego i trochę to smutne.