Basemap地图绘制_Python数据分析与可视化

Basemap是Matplotlib的一个子包，负责地图绘制。在数据可视化过程中，我们常需要将数据在地图上画出来。

比如说我们在地图上画出城市人口，飞机航线，军事基地，矿藏分布等等。这样的地理绘图有助于读者理解空间相关的信息。适用于有空间位置的数据集。

安装和使用

相对于其他工具Basemap用起来有点笨重，就算做点儿简单的可视化图也需要花费比预期更长的时间。

在处理比较复杂的地图可视化任务时，更现代的解决方案可能会更适用一些，比如leaflet 、Google Maps API。然而，Basemap 符合Python用户的使用习惯。

basemap并没有集成到matplotlib中，需要我们手动安装，basemap安装起来很简单。

安装命令：

pip install basemap

安装完并导入basemap工具箱后，只需要用几行代码就可以画出地理图：

import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap

plt.figure(figsize=(8, 8))
m = Basemap(projection='ortho', resolution=None, lat_0=50, lon_0=-100)
m.bluemarble(scale=0.5)

# 显式设置数据范围
plt.imshow(m.bluemarble(scale=0.5), origin='upper', vmin=0, vmax=1)

plt.show()

下面使用了 mill 投影方式，设置了地图的经纬度范围，绘制了海岸线、国家边界以及经纬度网格，也可以根据需要调整投影方式和绘制的内容

import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap

# 创建一个地理图
m = Basemap(projection='mill', llcrnrlat=-60, urcrnrlat=90, llcrnrlon=-180, urcrnrlon=180, resolution='c')

# 绘制海岸线和国家边界
m.drawcoastlines()
m.drawcountries()

# 绘制经纬度网格
m.drawparallels(range(-90, 91, 30), labels=[1, 0, 0, 0])
m.drawmeridians(range(-180, 181, 60), labels=[0, 0, 0, 1])

plt.title("World Map")
plt.show()

运用Basemap函数我们可以在绘图区域中绘制地理信息相关的图像，当参数 projection的值为'ortho'时，我们将得到一个如上图所示的地球仪截面。
将参数projection的值设置为lcc时，我们可以通过经纬度设置来得到某一区域的局部地图：

fig = plt.figure(figsize=(8, 8))
m = Basemap(projection='lcc', resolution=None,width=8E6,height=8E6,lat_0=45,lon_0=-100)m.etopo(scale=0.5, alpha=0.5)
# 将经纬度映射为 (x, y) 坐标，用于绘制图像
x, y = m(-122.3, 47.6)
plt.plot(x, y, 'ok', markersize=5)
plt.text(x, y, ' Seattle', fontsize=12)

这里使用了两个额外参数，它们用来表示地图中心的纬度（lat_0）和经度（ lon_0）。

地图投影

由于不可能把地表完美反映到二维平面上，所有的地图都是各种各样扭曲的产物，把这些扭曲的产物抹平到平面坐标系的过程，称为投影。

Basemap提供了几十种不同的投影的实现。

投影简写-全称对照：

下面我们对一常用的投影进行简单的演示。定义一个可以画带经纬线地图的简便方法：

def draw_map(m, scale=0.2):
    # 画地貌晕渲图
    m.shadedrelief(scale=scale)
    # 用字典表示经纬度
    lats = m.drawparallels(np.linspace(-90, 90, 13))
    lons = m.drawmeridians(np.linspace(-180, 180, 13))
    # 字典的键是plt.Line2D示例
    lat_lines = chain(*(tup[1][0] for tup in lats.items()))
    lon_lines = chain(*(tup[1][0] for tup in lons.items()))
    all_lines = chain(lat_lines, lon_lines)
    # 用循环将所有线设置成需要的样式
    for line in all_lines:
        line.set(linestyle='-', alpha=0.3, color='w')

圆柱投影是最简单的地图投影类型，纬度线与经度线分别映射成水平线与竖直线。
采用这种投影类型的话，赤道区域的显示效果非常好，但是南北极附近的区域就会严重变形。

fig = plt.figure(figsize=(8, 6), edgecolor='w')
m = Basemap(projection='cyl', resolution=None,
llcrnrlat=-90, urcrnrlat=90,
llcrnrlon=-180, urcrnrlon=180, )
draw_map(m)

这里basemap参数设置了左下角（llcrnr）和右上角（urcrnr）纬度（lat）和经度（lon）。不同的投影都有各种的优劣，大家之后可以多多尝试。

地图背景

basemap程序包中有许多实用的函数，可以画出各种地形的轮廓，如陆地、海洋、湖泊、河流、各国的政治分界线。

常用画图函数：

如果要使用边界特征，就必须设置分辨率。通过resolution来设置分辨率，取值为c（原始分辨率）、l（低分辨率）、i（中分辨率）、h（高分辨率）、f（全画质分辨率）。

来看看两种不同分辨率的绘制效果:

fig, ax = plt.subplots(1, 2, figsize=(12, 8))
for i, res in enumerate(['l', 'h']):
    m = Basemap(projection='gnom', lat_0=57.3, lon_0=-6.2,
        width=90000, height=120000, resolution=res, ax=ax[i])
    m.fillcontinents(color="#FFDDCC", lake_color='#DDEEFF')
    m.drawmapboundary(fill_color="#DDEEFF")
    m.drawcoastlines()
    ax[i].set_title("resolution='{0}'".format(res));
plt.show()

可以看出低分辨率不适合这个缩放，低分辨率适合呈现全局视角，而且加载速度比高分辨率更快。要呈现某一视角的适合，最好先从一个能快速呈现的分辨率开始，然后不断提高分辨率直到满意为止。

在地图上画数据

basemap还可以以地图为背景，在这上面画各种数据。basemap实例中许多方法都是与地图有关的函数。这些函数与标准matplotlib函数的用法类似，只是多了一个参数latlon。如果设置为true表示使用原来的经纬度坐标，不使用投影(x,y)坐标。

示例如下：

import pandas as pd
cities = pd.read_csv('california_cities.csv')
# 提取我们感兴趣的数据
lat = cities['latd'].values
lon = cities['longd'].values
population = cities['population_total'].values
area = cities['area_total_km2'].values
# 1. 绘制地图背景
fig = plt.figure(figsize=(8, 8))
m = Basemap(projection='lcc', resolution='h', 
            lat_0=37.5, lon_0=-119,
            width=1E6, height=1.2E6)
m.shadedrelief()
m.drawcoastlines(color='gray')
m.drawcountries(color='gray')
m.drawstates(color='gray')
# 2. 绘制城市数据的散点图，其中颜色反映人口
# 尺寸反映面积
m.scatter(lon, lat, latlon=True,
          c=np.log10(population), s=area,
          cmap='Reds', alpha=0.5)
# 3. 创建颜色条和图例
plt.colorbar(label=r'$\log_{10}({\rm population})$')
plt.clim(3, 7)
# 使用虚拟的点生成图例
for a in [100, 300, 500]:
    plt.scatter([], [], c='k', alpha=0.5, s=a,
                label=str(a) + ' km$^2$')
plt.legend(scatterpoints=1, frameon=False,
           labelspacing=1, loc='lower left');