简介:面向2024年华为杯研究生数学建模竞赛D题“大数据驱动的地理综合问题”,这套备赛资源围绕赛题各问展开,适合正在冲刺研赛、需要快速建立解题框架的队伍,也可作为大数据地理建模方向毕业设计的参考资料。内容覆盖问题拆解、思路分析、数学模型建立、算法实现到论文撰写全过程:每问均有超详细思路说明,模型部分给出完整建立过程,代码提供Matlab与Python双版本实现,便于在不同环境中运行与结果对比;配套的成品论文有助于梳理正文结构、图表与公式表达,可直接借鉴其论述节奏和章节安排。资源包以7z格式压缩,大小约10.27MB,内含原创思路、参考代码、PDF版思路文档等,信息密度较高,便于按需查阅和对照练习,尤其适合备赛时间紧张、需要快速入门的队伍。目前已有231人学习,若在处理数据、构建模型或解释结果时卡壳,可依照文件中从思路到代码的步骤逐段拆解,省去大量泛读资料的时间,把精力集中在D题核心模型的优化上。
1. 从数据到决策:D题地理综合问题的破题顺序
2024华为杯D题“大数据驱动的地理综合问题”看着像地理题,实际是典型的多源数据融合和空间建模。很多队伍第一天就开跑聚类,最后发现坐标系没对齐,聚合出来的网格特征全是错的。地理综合的核心不是高深算法,而是把人口、POI、遥感、气象等不同尺度的数据统一到同一空间单元,再通过加权、降维和聚类识别出有地理意义的综合分区。我按自己拆这道题的顺序写:数据清洗与网格对齐、特征工程与综合指数、空间聚类与调参、结果输出与竞赛技巧。每一步都给出可运行的Python求解代码,整套流程在普通笔记本上也能跑完。
2. 数据预处理与多源融合:从清洗到网格对齐
地理综合题的数据量通常不大,但脏数据不少。原始数据可能给你一张几十万行的csv,每行是一个网格单元,也可能给你几个shapefile和一个tif,要你自己做空间连接。第一步先做“体检”:用info()看类型,用describe()看量纲,用isnull().mean()看缺失比例。这一阶段的目标不是追求完美,而是让所有数据在同一个网格体系下对齐,后面所有分析才有意义。
2.1 多源数据清洗的核心流程
不同数据类型有不同问题,处理方法也不一样。我在处理时按下面这张表分工:
| 数据类型 | 常见问题 | 处理方式 |
|---|---|---|
| 站点观测数据 | 时间断档、短时异常 | 线性插值 + IQR截断 |
| POI数据 | 重复、坐标漂移 | 去重后按网格聚合 |
| 栅格数据 | 投影不一致、像元大小不一 | 重采样到统一网格 |
| 矢量路网 | 字段缺失、拓扑错误 | 计算密度后连接网格 |
地理数据的缺失值不能直接填均值,因为相邻网格往往有空间自相关。我一般用KNNImputer,用邻近网格的特征做插补。下面代码演示对三个连续型特征做邻域插补:
import pandas as pd from sklearn.impute import KNNImputer df = pd.read_csv('grid_features.csv') cols = ['nightlight', 'pop_density', 'road_density'] imputer = KNNImputer(n_neighbors=5, weights='distance') df[cols] = imputer.fit_transform(df[cols])n_neighbors=5表示参考最近5个网格,weights='distance'让距离近的样本权重更大。如果某个特征缺失率超过30%,建议直接剔除,插补出来的值会带偏后续PCA。实际跑题时,我曾经因为一个缺失率40%的“夜间灯光”特征补出一片带状噪声,后来发现是原始数据年份不一致,删掉反而更稳。
异常值用IQR截断,但要注意地理边界:一个网格的夜间灯光值可能因为行政区划差异很大,全局截断会把真实的高值区抹掉。我通常分区域做一次IQR,代码:
def clip_by_region(df, region_col='city', value_col='nightlight'): q1 = df.groupby(region_col)[value_col].transform('quantile', 0.25) q3 = df.groupby(region_col)[value_col].transform('quantile', 0.75) iqr = q3 - q1 lower = q1 - 1.5 * iqr upper = q3 + 1.5 * iqr return df[value_col].clip(lower, upper)transform会把分组统计量广播回每一行,避免手写循环。1.5倍IQR是比较常规的阈值,如果数据偏态严重,可以放宽到3倍。这个方法对人口密度、POI计数同样适用,只是阈值因子需要看直方图微调。
2.2 坐标系统一与网格聚合
多源数据要能合并,必须先统一坐标系。我用geopandas读矢量,统一转成EPSG:4326,但做距离计算时要再转投影坐标,否则以度为单位的结果会扭曲。栅格数据用rasterio读出来后,检查crs和transform,用reproject对齐网格。这里举一个POI聚合到网格的例子:
import geopandas as gpd grid = gpd.read_file('grid_1km.shp') poi = gpd.read_file('poi.shp') grid = grid.to_crs('EPSG:4326') poi = poi.to_crs('EPSG:4326') joined = gpd.sjoin(poi, grid, how='left', predicate='within') counts = joined.groupby('grid_id').size().reset_index(name='poi_count') grid = grid.merge(counts, on='grid_id', how='left') grid['poi_count'] = grid['poi_count'].fillna(0)sjoin的predicate='within'表示POI落在网格内部。如果POI数据有几十万行,空间连接会很慢,可以在连接前把经纬度转成float32减少内存。对于路网密度,可以用length / area,也可以像这样把路网线要素与网格相交后按长度累加。常见错误是直接计算POI点个数,忽略了网格面积差异,规范做法是先求出每个网格的POI密度再进模型。
在跨坐标系融合和网格聚合之后,马上做一次“数量核对”:把聚合后的网格数量和原始网格数核对,防止空间连接产生重复或丢失。地理综合题后面所有结果都建立在这个阶段,这一块出错,再好的模型都救不回来。
3. 特征工程与综合指数构建:让散点数据变成可分析的地理标签
数据对齐之后,常见的状态是每个网格有十多个特征:夜间灯光、人口密度、POI密度、路网密度、空气质量、气温、土地利用比例等。这些特征量纲不同,彼此高度相关。直接扔进聚类模型会让数值大的变量主导距离,需要先做特征筛选和标准化。
3.1 特征选择与归一化
先做相关性矩阵,把相关系数超过0.9的变量删掉一组。地理特征之间本来就强相关,比如夜间灯光和路网密度可能接近0.85,如果都保留,相当于给同一信息加了双倍权重。我用随机森林重要性辅助判断:
from sklearn.ensemble import RandomForestRegressor X = df[feature_cols].dropna() y = df['target_score'] # 可以是人口密度或其它引导变量 rf = RandomForestRegressor(n_estimators=200, random_state=42) rf.fit(X, y) imp = pd.Series(rf.feature_importances_, index=feature_cols).sort_values(ascending=False) print(imp)n_estimators=200是树的数量,随机森林对特征重要性估计比较稳定,不需要调太多。如果题目给的是分类标签,换成RandomForestClassifier。选前80%特征即可,别把重要性的台阶切得太细。
归一化我用MinMaxScaler,因为综合指数最终要映射到0到1区间,方便解释。StandardScaler适合后续聚类或PCA,但最终结果不如MinMax直观。代码:
from sklearn.preprocessing import MinMaxScaler scaler = MinMaxScaler() scaled = scaler.fit_transform(X[selected_cols]) scaled_df = pd.DataFrame(scaled, columns=selected_cols)fit_transform返回numpy数组,转成DataFrame是为了保留列名。注意fit只在训练集上调用,测试集或新数据用transform。
3.2 地理加权主成分分析
普通PCA忽略空间位置,但这道题叫“地理综合”,结果必须带空间结构。常见做法是把经纬度作为附加特征加入PCA,但仍然不是严格的地理加权。更稳妥的两步:先做主成分分析降维,再对主成分做空间自相关检验,确保提取出的综合因子不是随机噪点。
from sklearn.decomposition import PCA pca = PCA(n_components=0.9) pc = pca.fit_transform(scaled_df) print(pca.explained_variance_ratio_.cumsum())n_components=0.9表示保留累计90%方差的主成分。可以用碎石图判断,如果第3个主成分贡献率已经很低,可以只取前2个。地理综合问题通常有1到2个主导因子,比如“城市发展强度”和“生态本底”。
空间自相关检验可以用Moran's I,如果值接近0说明主成分没有空间结构,需要检查数据是否在网格上错位:
import esda from libpysal.weights import Queen w = Queen.from_dataframe(grid, use_index=True) mi = esda.moran.Moran(pc[:, 0], w) print(f"Moran's I = {mi.I:.3f}, p = {mi.p_sim:.3f}")Queen权重把共边网格视为邻居,适合网格数据。Moran's I大于0.3且p值小于0.05,说明有明显的空间聚类,地理综合因子的提取是有效的。我见过一些队伍直接从原始特征聚类,完全不检查空间结构,最后画出来的图像椒盐噪声一样碎,原因就在这里。
3.3 综合指数计算
地理综合指数最常用的是熵权法:根据指标离散程度确定权重,离散程度越大权重越高。相比主观AHP,它不需要专家打分,适合竞赛快速迭代。下面实现熵权法:
import numpy as np import pandas as pd def entropy_weight(x): x = np.asarray(x, dtype=float) # 平移避免0值 x = x - x.min(axis=0) + 1e-6 z = x / x.sum(axis=0)[None, :] k = 1 / np.log(z.shape[0]) e = -k * (z * np.log(z)).sum(axis=0) w = (1 - e) / (1 - e).sum() return w weights = entropy_weight(scaled_df.values) idx = pd.Series(weights, index=selected_cols).sort_values(ascending=False)代码中先做非负平移,避免log(0)。z是每个指标的概率分布,e是信息熵,权重等于信息冗余度的归一化。得到的weights可以直接和标准化值矩阵点乘得到综合得分:
score = scaled_df.values @ weights grid['geo_index'] = score最终geo_index就是每个网格的综合地理指数,范围在0到1之间,值越大表示综合发展水平或综合强度越高。
表格示例:
| 指标 | 信息熵 e | 权重 w |
|---|---|---|
| 夜间灯光 | 0.912 | 0.213 |
| 人口密度 | 0.934 | 0.158 |
| POI密度 | 0.887 | 0.271 |
| 路网密度 | 0.955 | 0.107 |
| 空气质量 | 0.979 | 0.051 |
实际权重会随数据变化,但规律通常是离散程度越高的指标权重越大。这一步做出来的综合指数可以直接用于后续功能区聚类和回归分析。
4. 模型求解与调参:从聚类到功能区识别
有了综合指数和标准化特征,下一步是把网格分成几个有明确地理意义的类别。地理综合题目里“分区分带”几乎是必考步骤。聚类算法选择很关键,不能直接拿K-Means默认参数一把梭。
4.1 空间聚类模型选择
常用三种算法对比如下:
| 算法 | 适合场景 | 缺点 |
|---|---|---|
| K-Means | 数据均匀、簇近似球形 | 对噪声敏感,需指定K |
| DBSCAN | 不规则簇、有噪声 | 密度相差大时效果差 |
| 层次聚类 | 需要嵌套结构 | 样本多时慢 |
地理网格数据通常没有明显噪声,K-Means最省心。但它不感知空间位置,所以我直接把归一化后的经纬度拼进特征矩阵,让空间邻近性参与聚类。注意经纬度要先投影成平面坐标(如EPSG:3857或UTM),否则纬度1度和经度1度距离不同,聚类结果会偏向纬度方向。
4.2 聚类实现与超参数选择
K的取值用轮廓系数在合理范围内扫一遍,不要直接拍脑袋定3类。代码:
from sklearn.cluster import KMeans from sklearn.metrics import silhouette_score X = np.hstack([scaled_df.values, projected_xy[['x', 'y']].values]) best_k, best_s = 2, -1 for k in range(2, 8): km = KMeans(n_clusters=k, n_init=10, random_state=42) labels = km.fit_predict(X) s = silhouette_score(X, labels) print(f"k={k}, silhouette={s:.4f}") if s > best_s: best_k, best_s = k, s final_labels = KMeans(n_clusters=best_k, n_init=10, random_state=42).fit_predict(X) grid['cluster'] = final_labelsn_init=10让K-Means用10个随机初值取最优,避免局部最优。轮廓系数的范围是-1到1,地理数据能到0.35以上已经算清晰。如果所有k的轮廓系数都低于0.25,先回到特征工程,把没有区分度的变量删掉。
4.3 空间自相关验证与稳定性检查
聚类结果应该满足“同区内相似、不同区间差异大”,而且最好具有空间连续性。用Moran's I验证:
mi = esda.moran.Moran(final_labels, w) print(f"cluster Moran's I = {mi.I:.3f}")另外要做稳定性测试:用不同random_state多跑几次,计算调整兰德指数(ARI)。如果两次结果差异很大,说明聚类不够稳定。
from sklearn.metrics import adjusted_rand_score lab1 = KMeans(n_clusters=best_k, n_init=10, random_state=0).fit_predict(X) lab2 = KMeans(n_clusters=best_k, n_init=10, random_state=1).fit_predict(X) print(adjusted_rand_score(lab1, lab2))ARI接近1表示几乎一致。如果低于0.8,可以考虑减少特征数量或改用层次聚类。聚类后的每个类别要做特征画像:计算各类的均值、中位数,判断哪个是市中心、哪个是工业区、哪个是生态区。这一步是论文“结果分析”部分的核心素材。
5. 结果输出与竞赛实用技巧
5.1 可视化与制图
用geopandas直接打印聚类专题图:
import matplotlib.pyplot as plt fig, ax = plt.subplots(figsize=(8, 8)) grid.plot(column='cluster', ax=ax, cmap='tab10', legend=True, edgecolor='white', linewidth=0.1) plt.axis('off') plt.savefig('cluster_map.png', dpi=300)dpi=300保证论文插图清晰。如果图例带数字,最好修改列名来显示功能区名称。综合指数用等间距分箱着色,而不是连续色阶,方便读者识别梯度。
5.2 论文表格生成技巧
Pandas可以输出LaTeX表格,竞赛时写论文很快:
grid.groupby('cluster')[['geo_index', 'nightlight', 'pop_density']].mean().round(3).to_latex('cluster_profile.tex')5.3 常见坑与检查清单
坐标系问题:sjoin前确保两个图层crs一致,否则会跑出空连接;缺失值问题:不要用全局均值填充,用邻域插值;网格面积问题:不同大小网格要加面积权重;PCA问题:标准化后贡献率才有意义;聚类问题:特征里加入经纬度前先投影。最后把grid_id重复检测和坐标系核对写成一个函数,每处理一个阶段就调用一次,确保merge后行数没有变化。这套流程我已经在多个城市大数据分析项目里验证过,拿来处理华为杯D题或者类似的大数据毕业设计题目,都能直接套用。
本文还有配套的精品资源,点击获取