当前位置：首页>python>使用 Python 绘制城市间专利合作申请数量网络图(一)

使用 Python 绘制城市间专利合作申请数量网络图(一)

2026-06-30 21:41:00

由于借助 AI 工具学习编程已经变得非常容易了，因此之后的课程就不再默认进行视频讲解了，如果特别需要视频讲解也可以联系李老师预约讲解～讲义材料学习过程中遇到的问题也可以及时与李老师联系。

前不久给大家分享了「1985～2024年省份间、城市间及区县间各类型专利合作数量面板数据」，在推文中我给大家展示了区县间、城市间和省份间的专利合作网络。今天我们以城市间专利申请合作网络为例讲解如何使用 Python 绘制这周空间网络图。

本次课程使用 Python 的 matplotlib + geopandas 复现了 R 语言中使用 ggplot2 基础图层绘制空间网络图的方法。

在学习本课程之前需要预先学习之前的课程：

使用 R 语言绘制历年中国省市区县地图（小地图版本+长版）: https://rstata.duanshu.com/#/brief/course/379d19770956478ebc4d919910fca74c

使用 R 语言绘制城市间专利合作申请数量网络图（一）: https://rstata.duanshu.com/#/brief/course/xxxxx

本课程中使用的地图数据也都是来自该课程。

使用 reticulate 创建与管理 Python 虚拟环境

在 R 中通过 reticulate 包来调用 Python，最好的实践是为项目创建一个专属的 Python 虚拟环境，将所需依赖隔离到独立空间，避免与系统 Python（如 Anaconda）发生版本冲突。

重要说明（避免"已初始化"报错）：reticulate 在 R 会话中只能绑定一次 Python——一旦某个 {python} 代码块运行，Python 解释器就被锁定，之后再调用 use_virtualenv() 会报错：
ERROR: The requested version of Python cannot be used, as another version has already been initialized.
因此，虚拟环境的激活必须在所有 {python} 代码块之前完成。本文档的解决方案是在 setup chunk 中通过 Sys.setenv(RETICULATE_PYTHON = ...) 提前锁定 Python 路径，这是 reticulate 选取 Python 的最高优先级入口。

安装 reticulate（仅首次）

# 设置 CRAN 镜像（knit 时 R 处于非交互模式，不会自动选择镜像）

options(repos = c(CRAN = "https://mirrors.tuna.tsinghua.edu.cn/CRAN/"))

# 仅在尚未安装时才安装，避免每次 knit 都重装

if (!requireNamespace("reticulate", quietly = TRUE)) {

install.packages("reticulate")

message("reticulate 安装完成！")

} else {

message("reticulate 已安装，版本：", packageVersion("reticulate"))

}

虚拟环境初始化原理（已在 setup chunk 中完成）

本文档的 setup chunk（隐藏运行）包含如下逻辑：

library(reticulate)

.venv_name <-".venv"

.venv_python <- virtualenv_python(.venv_name)

# 虚拟环境不存在时自动创建

if(!file.exists(.venv_python)){

virtualenv_create(.venv_name)

.venv_python <- virtualenv_python(.venv_name)

}

# 通过环境变量抢先锁定 Python（优先级最高，早于任何 {python} chunk）

Sys.setenv(RETICULATE_PYTHON = .venv_python)

use_virtualenv(.venv_name, required =TRUE)

这样做的关键在于：knitr 在处理第一个 {python} chunk 时，reticulate 已经通过 RETICULATE_PYTHON 环境变量知道要使用 .venv，不会再去碰 Anaconda。

数据读取与预处理

加载 Python 包

import pandas as pd

import numpy as np

import geopandas as gpd

import matplotlib.pyplot as plt

import matplotlib.font_manager as fm

import matplotlib.patches as mpatches

from matplotlib.lines import Line2D

from matplotlib.colors import to_rgba

import os

# 加载中文字体

font_path = "LXGWWenKai-Regular.ttf"

cnfont = fm.FontProperties(fname=font_path)

fm.fontManager.addfont(font_path)

font_name = fm.FontProperties(fname=font_path).get_name()

plt.rcParams['font.family'] = font_name

plt.rcParams['axes.unicode_minus'] = False

设置坐标系

中国地图通常使用 Albers 等面积投影：

MY_CRS = "+proj=aea +lat_0=0 +lon_0=105 +lat_1=25 +lat_2=47 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"

关键函数对照：

R	Python	说明
`sf::st_transform(mycrs)`	`gdf.to_crs(MY_CRS)`	坐标转换
`sf::st_point_on_surface()`	`gdf.geometry.representative_point()`	计算表面代表点
`sf::st_coordinates()`	`gdf.geometry.x / .y`	提取坐标
`sf::st_simplify(dTolerance=2000)`	`gdf.geometry.simplify(tolerance=2000)`	简化几何

读取专利合作数据

读取 2020 年城市间各类型专利合作数量统计数据：

df = pd.read_stata("2020年城市间各类型专利合作数量统计.dta")

数据清洗

重命名变量，并过滤掉三沙市（位置偏远且专利申请量极少）：

df = df.rename(columns={"城市1": "from", "城市2": "to",

"合作申请专利数量": "value"})

df = df[~df["from"].isin(["三沙市"]) & ~df["to"].isin(["三沙市"])]

df = df.reset_index(drop=True)

关键函数对照：

R (dplyr)	Python (pandas)	说明
`rename(from = 城市1, to = 城市2)`	`rename(columns={"城市1": "from", ...})`	重命名列
`filter(!from %in% "三沙市")`	`df[~df["from"].isin(["三沙市"])]`	条件筛选
`ungroup()`	无需操作	pandas 无 group 状态

统计每个城市的总合作量

# 每个城市的总合作量（等价于 R 的 group_by + summarise）

dfsum = df.groupby("from")["value"].sum().reset_index()

dfsum.columns = ["from", "value"]

dfsum

# 每条记录所在城市的总合作量（等价于 R 的 group_by + mutate(sum = sum(value))）

city_total = df.groupby("from")["value"].sum().reset_index()

city_total.columns = ["from", "sum"]

df = df.merge(city_total, on="from", how="left")

地图数据准备

读取地图数据

# 省界线（九段线、海岸线等）

provlinemap = gpd.read_file(

"chinaprov2021mini/chinaprov2021mini_line.shp", encoding="utf-8"

)

provlinemap = provlinemap[

~provlinemap["class"].str.contains("_", na=False) &

~provlinemap["class"].isin(["胡焕庸线", "秦岭-淮河线"])

]

provlinemap = provlinemap[["class", "geometry"]]

# 城市面

citymap = gpd.read_file(

"chinacity2021mini/chinacity2021mini.shp", encoding="utf-8"

)

citymap = citymap[citymap["省代码"].notna()].copy()

关键函数对照：

R (sf)	Python (geopandas)	说明
`read_sf("xxx.shp")`	`gpd.read_file("xxx.shp")`	读取矢量
`filter(!str_detect(class, "_"))`	`~df["class"].str.contains("_")`	正则筛选
`select(class)`	`df[["class", "geometry"]]`	选择列
`filter(!is.na(省代码))`	`df[df["省代码"].notna()]`	去除空值

计算城市质心

读取未经编辑的城市矢量数据，计算各城市的质心坐标（等价于 R 的 st_point_on_surface）：

此处代码需下载讲义材料查看～

匹配合作线起点和终点坐标

df = df.merge(city_centroiddf, left_on="from", right_on="市", how="left")

df = df.rename(columns={"X": "fromlon", "Y": "fromlat", "省": "prov"})

df = df.drop(columns=["市"])

df = df.merge(city_centroiddf, left_on="to", right_on="市", how="left")

df = df.rename(columns={"X": "tolon", "Y": "tolat"})

df = df.drop(columns=["市", "省"], errors="ignore")

简化地图数据

city_simplified = city.geometry.simplify(tolerance=2000, preserve_topology=True)

city = city.set_geometry(city_simplified)

准备散点和标签数据

由于城市数量众多，需要筛选出最重要的城市进行显示：

前 100 个城市绘制散点（按总合作量排序）
前 10 个城市添加文本标签（避免标签过于密集）

# 前100个城市用于散点（等价于 R 的 slice(1:100)）

first100df = (

df.groupby("from")

.agg(sum=("value", "sum"),

fromlon=("fromlon", "first"),

fromlat=("fromlat", "first"))

.reset_index()

.sort_values("sum", ascending=False)

.head(100)

.reset_index(drop=True)

)

# 前10个城市用于标签（等价于 R 的 slice(1:10)）

first10df = first100df.head(10).copy()

first10df

绘制网络图

配色方案

地图配色可以使用这个网站：https://tidyfriday.cn/colors

"小地图框格": "black",

"省份": "#b3b3b3"# R 的 gray30

}

line_widths = {

"九段线": 0.5,

"海岸线": 0.2,

"小地图框格": 0.2,

"省份": 0.3

}

# 前10个高亮城市的颜色

top10_colors = ["#fed439", "#709ae1", "#8a9197", "#d2af81", "#fd7446",

"#d5e4a2", "#197ec0", "#f05c3b", "#46732e", "#71d0f5"]

绘图核心代码

此处代码需下载讲义材料查看～

R vs Python 关键函数对照

本次课程涉及的 R → Python 核心函数映射：

R (ggplot2/sf)	Python (matplotlib/geopandas)	说明
`geom_sf(fill=NA, color="gray10")`	`gdf.plot(ax=ax, facecolor="none", edgecolor="#1a1a1a")`	绘制面边界
`geom_curve(curvature=0.3)`	`ax.annotate(..., connectionstyle="arc3,rad=0.3")`	弧线连接
`scale_linewidth_continuous(range=c(0.1, 0.5))`	`0.1 + 0.4 * (x - min) / (max - min)`	线宽映射
`scale_size_continuous(range=c(0.01, 7))`	`0.01 + 6.99 * (x - min) / (max - min)`	散点大小映射
`scale_alpha_manual(values=c(rep(1,10), rep(0.2, n-10)))`	`alpha = 1.0 if city in top10 else 0.2`	透明度控制
`ggrepel::geom_text_repel()`	`adjustText.adjust_text()`	避免重叠文本标签
`annotation_scale()`	手绘比例尺	比例尺
`annotation_north_arrow()`	`ax.annotate(..., arrowprops=dict(arrowstyle="->"))`	指北针
`guide_legend(direction="horizontal")`	`ax.legend(ncol=2)`	图例布局
`st_point_on_surface()`	`gdf.geometry.representative_point()`	表面代表点
`st_coordinates()`	`gdf.geometry.x / .y`	提取坐标

省份间和区县间网络图

附件中也提供了省份和区县版本的绘图代码（main2.py 和 main3.py），感兴趣的小伙伴也可以拿来参考。主要的区别在于：

数据源不同：省份版使用 2020年省份间各类型专利合作数量统计.dta，区县版使用 2020年区县间各类型专利合作数量统计.dta
ID 字段不同：省份版用 省1/省2 匹配省份名称，区县版用 县代码1/县代码2 匹配区县代码
地图底图不同：省份版使用 chinaprov2021mini，区县版使用 chinacounty2021mini
线宽范围不同：省份版使用 range=c(0.2, 1.5)，区县版使用 range=c(0.1, 0.5)
透明度不同：区县版非 top10 城市使用 alpha=0.05（因为区县数量更多，需要更低的透明度）
标签内容不同：区县版标签使用 市+县 格式，而非单纯的区县名称

如何参加课程？

是不是感觉很硬核！欢迎报名 RStata 培训班获取全部课程和以会员价获取数据资料（10元/份）详情可阅读这篇推文：数据处理、图表绘制、效率分析与计量经济学如何学习～

详情可点击阅读原文进入 RStata 学院了解（从首页的会员卡专区即可查看和购买会员卡）。

更多关于 RStata 培训班的信息可添加微信号 r_stata2 咨询：

附件下载（点击文末的阅读原文即可跳转）：
https://rstata.duanshu.com/#/brief/course/23c96f80df094f64a6e6fbafd03b3043

本文来自网友投稿或网络内容，如有侵犯您的权益请联系我们删除，联系邮箱：wyl860211@qq.com 。

绘图核心代码

R vs Python 关键函数对照

省份间和区县间网络图

如何参加课程？

使用 Python 绘制城市间专利合作申请数量网络图(一)