はじめに
この実験では、ハーベルシン距離メトリック（つまり、緯度/経度の点間の距離）に基づいて構築されたボールツリーを使用して、地理空間データに対する近傍ベースの照会（特に核密度推定）の例を示します。データセットは Phillips et. al. (2006) によって提供されています。この例では、ベースマップライブラリを使用して南米の海岸線と国境を描画します。
VM のヒント
VM の起動が完了したら、左上隅をクリックしてノートブックタブに切り替え、Jupyter Notebook を使って練習しましょう。
時々、Jupyter Notebook が読み込み終了するまで数秒待つ必要がある場合があります。Jupyter Notebook の制限により、操作の検証を自動化することはできません。
学習中に問題に遭遇した場合は、Labby にお問い合わせください。セッション後にフィードバックを提供してください。私たちは迅速に問題を解決いたします。
必要なライブラリをインポートする
この実験で必要なライブラリをインポートするのが最初のステップです。この実験では、numpy、matplotlib、fetch_species_distributions、および KernelDensity ライブラリを使用します。
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import fetch_species_distributions
from sklearn.neighbors import KernelDensity
データを読み込む
次のステップは、Phillips et. al. (2006) によって提供されたデータを読み込むことです。このデータセットには、核密度推定を示すために使用する 2 つの種が含まれています。
data = fetch_species_distributions()
species_names = ["Bradypus Variegatus", "Microryzomys Minutus"]
データを準備する
ここでは、核密度推定のためのデータを準備します。データセットから緯度と経度の情報を抽出し、ラジアンに変換します。
Xtrain = np.vstack([data["train"]["dd lat"], data["train"]["dd long"]]).T
ytrain = np.array(
[d.decode("ascii").startswith("micro") for d in data["train"]["species"]],
dtype="int",
)
Xtrain *= np.pi / 180.0 ## Convert lat/long to radians
グリッドを構築する
ここでは、バッチオブジェクトから地図グリッドを構築します。これを達成するためにconstruct_grids関数を使用します。
def construct_grids(batch):
"""Construct the map grid from the batch object
Parameters
----------
batch : Batch object
The object returned by :func:`fetch_species_distributions`
Returns
-------
(xgrid, ygrid) : 1-D arrays
The grid corresponding to the values in batch.coverages
"""
## x,y coordinates for corner cells
xmin = batch.x_left_lower_corner + batch.grid_size
xmax = xmin + (batch.Nx * batch.grid_size)
ymin = batch.y_left_lower_corner + batch.grid_size
ymax = ymin + (batch.Ny * batch.grid_size)
## x coordinates of the grid cells
xgrid = np.arange(xmin, xmax, batch.grid_size)
## y coordinates of the grid cells
ygrid = np.arange(ymin, ymax, batch.grid_size)
return (xgrid, ygrid)
## Call the function and store the results in xgrid and ygrid
xgrid, ygrid = construct_grids(data)
データグリッドを準備する
コントーアプロット用のデータグリッドを設定します。これを達成するためにconstruct_grids関数を使用します。
X, Y = np.meshgrid(xgrid[::5], ygrid[::5][::-1])
land_reference = data.coverages[6][::5, ::5]
land_mask = (land_reference > -9999).ravel()
xy = np.vstack([Y.ravel(), X.ravel()]).T
xy = xy[land_mask]
xy *= np.pi / 180.0
南アメリカの地図を描画する
ここでは、各種の分布を伴う南アメリカの地図を描画します。
fig = plt.figure()
fig.subplots_adjust(left=0.05, right=0.95, wspace=0.05)
for i in range(2):
plt.subplot(1, 2, i + 1)
print(" - computing KDE in spherical coordinates")
kde = KernelDensity(
bandwidth=0.04, metric="haversine", kernel="gaussian", algorithm="ball_tree"
)
kde.fit(Xtrain[ytrain == i])
Z = np.full(land_mask.shape[0], -9999, dtype="int")
Z[land_mask] = np.exp(kde.score_samples(xy))
Z = Z.reshape(X.shape)
levels = np.linspace(0, Z.max(), 25)
plt.contourf(X, Y, Z, levels=levels, cmap=plt.cm.Reds)
if basemap:
print(" - plot coastlines using basemap")
m = Basemap(
projection="cyl",
llcrnrlat=Y.min(),
urcrnrlat=Y.max(),
llcrnrlon=X.min(),
urcrnrlon=X.max(),
resolution="c",
)
m.drawcoastlines()
m.drawcountries()
else:
print(" - plot coastlines from coverage")
plt.contour(
X, Y, land_reference, levels=[-9998], colors="k", linestyles="solid"
)
plt.xticks([])
plt.yticks([])
plt.title(species_names[i])
plt.show()
まとめ
この実験では、地理空間データに対する核密度推定を行う方法を学びました。この手法を示すために、Phillips et. al. (2006) のデータセットを使用しました。また、basemap ライブラリを使って、各種の分布を伴う南アメリカの地図を描画する方法も学びました。