[PYTHON] Visualisez et comprenez le maillage régional du Japon sur une carte

Qu'est-ce qu'un maillage régional?

** Le maillage régional ** est un maillage qui a presque la même forme que la région en fonction de la latitude et de la longitude. Un maillage est approximativement carré et est défini dans différentes tailles. Un bref résumé des choses que vous utilisez souvent

Référence: Bureau des statistiques, ministère de l'Intérieur et des Communications | Statistiques du maillage régional

Le système de code est fixe pour le maillage, et il est possible de calculer à partir de la latitude et de la longitude (Référence: Système de maillage standard et code) .. Cependant, la définition de chaque code de maillage est également fournie sous forme de données (e-Stat | Statistics on Map (Statistical GIS). -search? page = 1 & type = 2)). Dans cet article, je voudrais comprendre en utilisant ceci pour visualiser le maillage régional.

Bibliothèque de logiciels

Python (distribution Anaconda / miniconda) est utilisé. Les bibliothèques requises peuvent être installées à partir de la commande suivante.

conda install -c conda-forge jupyter geopandas descartes shapely \
  matplotlib pip requests pillow chardet mplleaflet && \
  pip install tilemapbase

Les bibliothèques suivantes sont caractéristiques de cet objectif.

Obtenez des données

Cette fois, nous utiliserons les données du maillage primaire 5339. Le fichier peut être téléchargé directement depuis ce lien ou Téléchargez le fichier de forme "M5339" à partir de la page e-Stat. Je vais. Puisqu'il s'agit d'un fichier Zip, développez-le et s'il a la structure suivante, vous êtes prêt à commencer.

└── QDDSWQ5339
    ├── MESH05339.dbf
    ├── MESH05339.prj
    ├── MESH05339.shp
    └── MESH05339.shx

Lire les données

Le chargement des fichiers de forme est très facile avec la bibliothèque geopandas.

import geopandas as gpd

x = gpd.read_file("QDDSWQ5339/MESH05339.shp")
print(x.shape)
x.head()

(100800, 8)
     KEY_CODE MESH1_ID MESH2_ID MESH3_ID MESH4_ID MESH5_ID  OBJ_ID  \
0  5339000011     5339       00       00        1        1       1   
1  5339000012     5339       00       00        1        2       2   
2  5339000013     5339       00       00        1        3       3   
3  5339000014     5339       00       00        1        4       4   
4  5339000021     5339       00       00        2        1       5   

                                            geometry  
0  POLYGON ((139.00312 35.33333, 139.00000 35.333...  
1  POLYGON ((139.00625 35.33333, 139.00312 35.333...  
2  POLYGON ((139.00312 35.33542, 139.00000 35.335...  
3  POLYGON ((139.00625 35.33542, 139.00312 35.335...  
4  POLYGON ((139.00937 35.33333, 139.00625 35.333...

Agrégation de polygones

Puisque les données cette fois sont les données du maillage du 5ème ordre, par exemple, l'ensemble du maillage du 1er ordre "5339" n'est défini nulle part. Pour un maillage grossier, il peut être obtenu en agrégeant des polygones de 5ème ordre. La méthode dissolve de geopandas.GeoDataFrame est utile pour cela.

La fonction suivante agrège les données du maillage de 5ème ordre aux données d'ordre inférieur.

def aggregate_mesh(x, level):
    tmp = x.copy()
    code_len = [4, 6, 8, 9, 10][level-1]
    tmp["key"] = tmp["KEY_CODE"].str.slice(0, code_len)
    tmp = tmp[["key", "geometry"]]
    tmp = tmp.dissolve(by="key")
    return tmp

--code_len: longueur du code avec la granularité spécifiée --Spécifiez le code de granularité spécifié dans tmp [" key "] --Agréger les polygones avec .dissolve --. Dissoudre prend une somme pour les informations géographiques. Je ne l'ai pas utilisé cette fois, mais s'il y a d'autres colonnes, il fera des calculs comme groupby.aggregate en même temps.

# aggregate to mesh level 1
mesh1 = aggregate_mesh(x, 1)
mesh1

      geometry
key                                                    
5339  POLYGON ((139.29062 35.33333, 139.28750 35.333...

C'est le résultat de l'agrégation dans le maillage principal. Cette fois, seul le maillage principal de "5339" est utilisé, il sera donc agrégé en une seule ligne.

Visualisation des informations géographiques du maillage (carte interactive)

Le module geopandas.plotting fournit des fonctions utiles pour visualiser les informations géographiques. Surtout pour les polygones, plot_polygon_collection est utile. Littéralement, c'est une fonction qui visualise plusieurs polygones.

import matplotlib.pyplot as plt

fig, ax = plt.subplots()
plot_polygon_collection(ax, mesh1.geometry, edgecolor="#121212", facecolor="#333333")

mesh1.png

Le résultat sera un graphique comme celui ci-dessus. Je ne sais pas o c'est. Utilisez la bibliothèque mplleaflet pour tracer ceci sur une carte. Cette bibliothèque offre la possibilité de dessiner des graphiques matplotlib sur une carte interactive.

import mplleaflet

mplleaflet.display(fig)  #Pour jupyter
# mplleaflet.show(fig)   #Autre

leaflet1.png

Lorsqu'elle est exécutée avec Jupyter, la carte sera affichée dans la sortie comme indiqué ci-dessus. Il peut être agrandi et déplacé comme une application cartographique normale. Si vous utilisez .show, une autre page du navigateur s'ouvrira et une carte similaire y sera affichée.

Vous pouvez voir que le maillage "5339" est un maillage qui comprend une partie de la région sud du Kanto depuis Tokyo. La raison pour laquelle ce n'est pas exactement carré est probablement parce que le maillage du 5ème ordre de la partie mer n'est pas défini.

De même, affichons le maillage secondaire.

mesh2 = aggregate_mesh(x, 2)

fig, ax = plt.subplots()
plot_polygon_collection(ax, mesh2.geometry, edgecolor="#121212", facecolor="#333333")
mplleaflet.display(fig)

leaflet2.png

Visualisation des informations géographiques du maillage (carte statique)

Bien que mplleaflet soit une fonctionnalité très utile, il présente également quelques faiblesses.

«Plus vous tracez, plus cela devient lourd.

Concernant le premier, il est difficile de dessiner, par exemple, si vous donnez des milliers à 10 000 diagrammes de dispersion. Comme alternative, essayez de dessiner sur une carte statique. Il semble y avoir quelques options, mais nous utiliserons la bibliothèque TileMapBase. Ce package prend une carte OpenStreetMap comme image et la dessine en arrière-plan. Vous pouvez utiliser librement les fonctions de matplotlib en superposant le dessin dessus. La faiblesse est que la syntaxe n'est pas aussi simple que «mplleaflet», et la conversion des coordonnées doit être effectuée par l'utilisateur.

Ce sera un peu plus long, mais c'est une fonction qui dessine un polygone sur la carte et dessine une chaîne de caractères au centre.

import tilemapbase
tilemapbase.init(create=True)

from shapely import geometry

def plot_on_map_with_text(geoseries, texts):
    # find the graph range
    rect = geoseries.total_bounds
    edgex = rect[2] - rect[0]
    edgey = rect[3] - rect[1]
    extent = tilemapbase.Extent.from_lonlat(
        rect[0]-edgex*0.3, rect[2]+edgex*0.3,
        rect[1]-edgey*0.3, rect[3]+edgey*0.3)
    extent = extent.to_aspect(1.0)

    fig, ax = plt.subplots(figsize=(8, 8), dpi=100)
    ax.xaxis.set_visible(False)
    ax.yaxis.set_visible(False)
    
    t = tilemapbase.tiles.build_OSM()
    plotter = tilemapbase.Plotter(extent, t, width=600)
    plotter.plot(ax, t)

    polygons = []
    centers = []
    bounds = geoseries.bounds
    for i in range(len(bounds)):
        # convert to the plottable scale
        minx, miny = tilemapbase.project(bounds['minx'][i], bounds['miny'][i])
        maxx, maxy = tilemapbase.project(bounds['maxx'][i], bounds['maxy'][i])
        polygons.append(
            geometry.box(minx, miny, maxx, maxy))
        centers.append([(minx + maxx)/2.0, (miny + maxy)/2.0])

    polygons = gpd.GeoSeries(polygons)
    plot_polygon_collection(ax, polygons,
                            edgecolor='#121212', facecolor='#000000', alpha=0.4)
    for center, txt in zip(centers, texts):
        ax.text(center[0], center[1], txt, fontdict={'color':'lightblue'})
    return fig, ax

--tilemapbase.init doit être exécuté une fois avant utilisation et crée un cache de l'image de la carte acquise pour éviter de récupérer les mêmes informations à plusieurs reprises. Cela semble être une considération pour ne pas déranger OpenStreetMap en répétant l'accès. --shapely est une bibliothèque qui fournit les fonctions des éléments de base de l'information géographique tels que les polygones, et geopandas en dépend également. Ici, il est utilisé pour redéfinir le polygone dont le système de coordonnées a été converti.

Tracez le maillage quadratique avec les deux derniers chiffres du code.

plot_on_map_with_text(mesh2.geometry, mesh2.index.str.slice(4, 6))

static2.png

En regardant cela, vous pouvez voir les règles du code de maillage. Avec l'extrémité sud-ouest comme «(0, 0)», le premier nombre augmente en allant vers le nord, et le deuxième nombre augmente en allant vers l'est. Le maillage secondaire divise le maillage principal en huit, de sorte que les nombres de 0 à 7 apparaissent.

Vous pouvez dessiner de la même manière après le 3ème ordre. Cependant, le graphique sera trop détaillé, alors limitez-le à une zone spécifique.

mesh3 = aggregate_mesh(x[x.MESH2_ID == "45"], 3)
plot_on_map_with_text(mesh3.geometry, mesh3.index.str.slice(6, 8))

static3.png

Les règles du code n'ont pas changé, et cette fois il est divisé en 10.

mesh4 = aggregate_mesh(x[(x.MESH2_ID == "45") & (x.MESH3_ID == "09")], 4)
plot_on_map_with_text(mesh4.geometry, mesh4.index.str.slice(8, 9))

static4.png Le maillage quaternaire divise le maillage quaternaire en deux parties verticalement et horizontalement, et attribue respectivement 1, 2, 3 et 4 codes au sud-ouest, sud-est, nord-ouest et nord-est.

mesh5 = aggregate_mesh(x[(x.MESH2_ID == "45") & (x.MESH3_ID == "09")], 5)
plot_on_map_with_text(mesh5.geometry, mesh5.index.str.slice(8, 10))

static5.png Les quatre maillages quaternaires sont divisés en cinq et affichés. Les règles de codage pour la 5e maille sont les mêmes que pour la 4e maille.

Recommended Posts

Visualisez et comprenez le maillage régional du Japon sur une carte
Visualisez les buts et les buts de la Ligue B sur une carte de chaleur
Folium: Visualisez les données sur une carte avec Python
Visualiser grib2 sur une carte avec python (matplotlib)
Visualisez les préfectures avec de nombreux itinéraires par préfecture sur une carte du Japon
Un mémo contenant Python2.7 et Python3 dans CentOS
Carte des informations de location sur une carte avec python
Comprenez les listes Python, les dictionnaires, etc.
Créer une carte Web en utilisant Python et GDAL
[Rails] Comment afficher plusieurs marqueurs sur Google Map et afficher une bulle lorsque vous cliquez dessus