Python: как создать карту choropleth из шейп файла Канады?
Моя цель здесь - создать choropleth map Канады в Python. Предположим, у меня есть словарь со значениями, относящимися к каждой канадской провинции/территории:
myvalues={'Alberta': 1.0,
'British Columbia': 2.0,
'Manitoba': 3.0,
'New Brunswick': 4.0,
'Newfoundland and Labrador': 5.0,
'Northwest Territories': 6.0,
'Nova Scotia': 7.0,
'Nunavut': 8.0,
'Ontario': 9.0,
'Prince Edward Island': 10.0,
'Quebec': 11.0,
'Saskatchewan': 12.0,
'Yukon': 13.0}
Теперь я хочу покрасить каждую провинцию на основе соответствующего значения в myvalues
, используя непрерывную цветовую палитру (например, оттенки красного). Как это сделать?
До сих пор я мог только строить канадские провинции/территорию в пределах matplotlib, но их формы появляются в уникальном цвете, и я не знаю, как изменить это в соответствии с числами в myvalues
(может быть, я нужно играть с patches
, но я не знаю, как).
Здесь вы можете найти шейп файл: http://www.filedropper.com/canadm1_1
И это мой код на сегодняшний день:
import shapefile
import matplotlib.pyplot as plt
import matplotlib.patches as patches
from matplotlib.patches import Polygon
from matplotlib.collections import PatchCollection
# -- input --
sf = shapefile.Reader("myfolder\CAN_adm1.shp")
recs = sf.records()
shapes = sf.shapes()
Nshp = len(shapes)
cns = []
for nshp in xrange(Nshp):
cns.append(recs[nshp][1])
cns = array(cns)
cm = get_cmap('Dark2')
cccol = cm(1.*arange(Nshp)/Nshp)
# -- plot --
fig = plt.figure()
ax = fig.add_subplot(111)
for nshp in xrange(Nshp):
ptchs = []
pts = array(shapes[nshp].points)
prt = shapes[nshp].parts
par = list(prt) + [pts.shape[0]]
for pij in xrange(len(prt)):
ptchs.append(Polygon(pts[par[pij]:par[pij+1]]))
ax.add_collection(PatchCollection(ptchs,facecolor=None,edgecolor='k', linewidths=.5))
ax.set_xlim(-160,-40)
ax.set_ylim(40,90)
Это изображение, которое я получаю до сих пор:
![введите описание изображения здесь]()
ИЗМЕНИТЬ
Я получаю решение должно быть в следующих строках:
cm = get_cmap('OrRd')
cccol = cm(1.*arange(Nshp)/Nshp)
Вышеупомянутый script создает массив cccol
, который на самом деле имеет эту форму:
array([[ 1. , 0.96862745, 0.9254902 , 1. ],
[ 0.99766244, 0.93356402, 0.84133796, 1. ],
[ 0.99520185, 0.89227221, 0.74749713, 1. ],
[ 0.99274125, 0.84306037, 0.64415227, 1. ],
[ 0.99215686, 0.78754327, 0.5740254 , 1. ],
[ 0.99186467, 0.71989237, 0.50508269, 1. ],
[ 0.98940408, 0.60670514, 0.39927722, 1. ],
[ 0.97304114, 0.50618995, 0.32915034, 1. ],
[ 0.94105344, 0.40776625, 0.28732027, 1. ],
[ 0.88521339, 0.28115341, 0.19344868, 1. ],
[ 0.8220992 , 0.16018455, 0.10345252, 1. ],
[ 0.73351789, 0.04207613, 0.02717416, 1. ],
[ 0.61959248, 0. , 0. , 1. ]])
Я не знаю, почему у него 4 столбца, но я полагаю, что если я каким-то образом свяжу эти значения с параметрами, указанными в values
dict, я могу решить эту проблему. Любые идеи?
РЕДАКТИРОВАТЬ 2
Я понял, что "трюк" находится в cccol = cm()
. Чтобы связать это с провинциями, я попытался назначить
cccol = cm(myvalues.values(i) for i in myvalues.keys())
чтобы (по моему мнению, по крайней мере) каждый цвет назначается на основе соответствующего ключа, и нет никаких ошибок. Проблема в том, что я получаю сообщение об ошибке:
TypeError: Cannot cast array data from dtype('O') to dtype('int32') according to the rule 'safe'
.
Как обойти это?
Ответы
Ответ 1
Это напрямую не отвечает на ваш вопрос, но, надеюсь, решает вашу проблему точно так же. Вы посмотрели GeoPandas? Он предоставляет простой API для работы с формами шейпов и их построения. Вы можете копировать свой код, включая прорисовку choropleth, всего несколькими строчками:
import geopandas as gpd
canada = gpd.read_file('CAN_adm1.shp')
canada.plot('myvalues', cmap='OrRd')
В этом примере предполагается, что ваш шейп файл имеет атрибут в каждой провинции, который содержит значения, которые вы хотите построить, и атрибут называется "myvalues". Если значения не сохранены в шейп файле, вы можете использовать canada.merge
для объединения вашей карты values
в GeoDataframe.
Одно предостережение: в это время GeoPandas не имеет простого способа построить легенду для цветов choropleth. (вопрос, указанный здесь)
Ответ 2
Запрос: переименуйте словарь values
в другое. Это имя заставило писать этот ответ гораздо сложнее.:)
Не тестировали это, но попробуйте:
color_numbers = values.values()
# assumes the provinces are listed in the same order in values as
# they are in the shape file
for nshp in xrange(Nshp):
ptchs = []
# ... code omitted ...
the_facecolor = [(color_numbers[nshp]-1)/(Nshp-1), 0, 0]; #1..13 -> 0..1, then add G=B=0.
# change the computation if the values in the values dictionary are no longer 1..13
ax.add_collection(PatchCollection(ptchs, facecolor=the_facecolor, edgecolor='k', linewidths=.5))
Получаемый вами результат имеет все синие патчи или [0,0,1]
. Поскольку эта строка не находится в cccol
, я не думаю, что проблема cccol
. Кроме того, добавленный код никогда не ссылался на cccol
после его создания! (Пожалуйста, добавьте ссылку на образец кода, с которого вы начали!:))
В любом случае, установка facecolor
должна помочь, насколько я знаю. Преобразуя запись values
в диапазон 0..1, затем делая цветные записи [R, G, B], должны давать вам оттенки красного цвета.
Ответ 3
Вы упомянули путаницу о cccol
, являющемся списком списков. Это список кортежей RGBA (красная, зеленая, синяя, альфа-прозрачность). Они представляют 13 "равноотстоящих" цветов от оранжевого до красного.
В вашем случае вам не нужны одинаково разнесенные цвета, а цвета, соответствующие myvalues
. Сделайте это:
cmap = matplotlib.cm.get_cmap('OrRd')
norm = matplotlib.colors.Normalize(min(myvalues.values()), max(myvalues.values()))
color_producer = matplotlib.cm.ScalarMappable(norm=norm, cmap=cmap)
Теперь color_producer
имеет метод to_rgba
, который принимает значения из myvalues
и преобразует их в правильные цвета. Normalize
устанавливает минимальный и максимальный диапазон myvalues
в крайние цвета красно-оранжевой цветовой карты.
Теперь, когда вы создаете каждую провинцию PatchCollection
, вы можете установить ее facecolor
в кортеж RGBA, возвращенный color_producer
:
# Change the province name passed as you iterate through provinces.
rgba = color_producer.to_rgba(myvalues['Manitoba'])
PatchCollection(ptchs, facecolor=rgba, edgecolor='k', linewidths=.5)