简介针对2024年全国大学生数学建模竞赛C题这份资源围绕农作物种植策略提供了一整套可落地的分析方案涵盖建模思路、程序源码和优秀论文适合参赛学生、赛事指导老师以及想入门数学建模的读者使用。代码按赛题问题拆分前两问均有独立求解脚本第三问综合运用线性回归、非线性回归、Pearson相关性检验和KMeans聚类算法对预期销售量、种植成本和销售单价进行统计建模并对地块—大棚做分类处理还附带整理好的地块信息、作物信息、相关性信息及2023年汇总补充数据基本还原了从数据处理到结果输出的完整链路。资源共28个文件包括12个Excel数据文件、8个Python程序脚本、2个PDF赛题与论文、2个Markdown说明文档及4个代码备份文件压缩包整体仅2.21MB轻量易用方便按需查阅。当前已有133人学习适合需要复现国赛获奖思路、参考建模代码或快速了解农作物种植优化问题的读者。1. 从赛题到代码一份完整 C 题方案的解构2024 年国赛 C 题把「农作物种植策略」做成了一个相当实际的优化问题村里有平旱地、梯田、山坡地、水浇地还有大棚十几种作物要在 2024—2030 年间逐年排产产量、售价、成本还会波动。这份资源里真正值钱的不是那篇论文而是六个按问题顺序编号的 Python 求解文件从 Q1_1.py、Q1_2.py 的线性规划起步到 Q2.py 的随机优化再到问题三的 pearson.py、线性/非线性回归脚本和 KMeans 聚类完整覆盖了确定性优化、不确定性建模、相关性检验和聚类四类数学建模模型。它适合两类人正在备战数学建模竞赛、想找一套可复现代码路径的选手以及用 Python 处理生产调度、想把优化工具链串起来落地的开发者。下文按赛题推进顺序逐层拆这套代码。2. 地块与作物的约束系统Q1_1.py 和 Q1_2.py 的线性规划实现2.1 从赛题数据到决策变量问题一的核心假设是未来产量、售价、成本相对稳定因此方案可以在同一套约束下逐年复用。数据基础是资源包里的「地块信息.xlsx」和「作物信息.xlsx」前者记录每块地的类型与面积后者记录每种作物的亩产量、种植成本、销售单价。建模的第一步是把自然语言描述转成线性规划的标准形式。决策变量定义为# x[i * n_crop j] 表示第 i 块地上种植第 j 种作物的面积 # 索引展开方式地块优先同一地块内的作物连续排列 n_land len(land) n_crop len(crop) x [0] * (n_land * n_crop)变量采用「地块优先 块内作物连续」的扁平排列是因为scipy.optimize.linprog只接受一维决策变量后续构造约束矩阵时需要按i * n_crop j这个索引规则定位。如果换成二维变量矩阵构造会直观一些但还原成线性规划标准形式时需要多一层 reshape反而容易出错。2.2 目标函数与约束矩阵的构造目标函数是总利润最大化。单块地上某作物的利润按「面积 × 亩产量 × 销售单价 − 面积 × 种植成本」计算因此目标系数向量里的每一项都是单位面积净利润。linprog只做最小化所以把利润向量整体取负。import pandas as pd from scipy.optimize import linprog def build_q1_model(land, crop): n_land, n_crop len(land), len(crop) c [] for _, ld in land.iterrows(): for _, cr in crop.iterrows(): profit cr[亩产量] * cr[销售单价] - cr[种植成本] c.append(-profit) # 最小化负利润 最大化利润 A_ub, b_ub [], [] # 约束每块地上所有作物的种植面积之和 该地块面积 for i in range(n_land): row [0] * (n_land * n_crop) for j in range(n_crop): row[i * n_crop j] 1 A_ub.append(row) b_ub.append(land.loc[i, 面积]) # 变量下界为 0上界为该地块面积保证单个变量不会超过地块总面积 bounds [(0, land.loc[i, 面积]) for i in range(n_land) for _ in range(n_crop)] res linprog(c, A_ubA_ub, b_ubb_ub, boundsbounds, methodhighs) return res这里的methodhighs是 SciPy 1.7 之后默认启用的 HiGHS 求解器适合中等规模的线性规划如果直接用旧代码里的methodsimplex在数千个变量的场景下迭代步数和数值误差都会明显变大。约束矩阵 A_ub 的每一行只在一块地对应的 n_crop 个位置上填 1本质是块内加和约束这在实际求解中比单独对每块地每种作物分别设限更接近赛题的「可以少种但不可以多种」逻辑。2.3 Q1_2.py 在 Q1_1 之上叠加的时间维约束Q1_1.py 求解的是稳定假设下的一期最优布局而 Q1_2.py 解决的是问题一的第二问方案要逐年滚动到 2030 年同一块地上不能连续多年种植同一种作物。此时决策变量需要扩成三维x[t][i][j]表示第 t 年第 i 块地种第 j 种作物的面积并额外追加重茬约束。T 7 # 2024 到 2030 def add_rotation_constraint(A_ub, b_ub, n_land, n_crop, T, forbid_pairs): for t in range(1, T): # 从第二年开始检查上一年的种植情况 for i in range(n_land): for j in forbid_pairs.get(i, []): row [0] * (T * n_land * n_crop) # 上一年的面积 这一年的面积 上一年面积的阈值 row[(t-1) * n_land * n_crop i * n_crop j] 1 row[t * n_land * n_crop i * n_crop j] 0.2 A_ub.append(row) b_ub.append(0) return A_ub, b_ub这段约束的含义是如果第 t−1 年某地块种了易重茬障碍的作物第 t 年的同作物面积不能超过上一年的 0.2 倍相当于允许少量留茬但不允许连续性大规模种植。forbid_pairs需要根据作物信息里的「重茬敏感」标记手动配置不同作物对连作的容忍度差异很大这属于赛题数据里不会直接给出、但实际求解时必须补全的工程细节。Q1_1 和 Q1_2 的关系可以简单理解成前者是「一年定终身」后者是「七年滚动排产」。从代码量看Q1_2 只多了时间索引和轮作约束但求解规模扩大了近 7 倍这也是我建议用highs而不是老求解器的直接原因。3. 不确定性建模Q2.py 与改进粒子群优化算法IPSO3.1 为什么线性规划在问题二失效问题二把产量、售价、成本同时变成随机变量。如果继续用线性规划就得把每个随机参数替换成期望值但期望值解在真实波动下往往违背面积约束或产生明显亏损。原因是线性规划天然假设参数确定一旦目标函数里的系数带随机性最优解对参数扰动非常敏感尤其是产量和价格呈负相关时丰年产量高但价格低直接用均值建模会系统性高估收益。资源包里的论文标题给出了作者团队的解法基于粒子群优化算法的华北农作物种植策略模型。代码里对应的是 Q2.py它引入三个情景——丰年、平年、歉年每个情景有独立的产量折减系数和价格波动系数目标函数改为多年期望利润。import numpy as np from sko.PSO import PSO def expected_profit(x, land, crop, years7): X x.reshape(years, len(land), len(crop)) scenarios {丰年: (1.10, 0.95), # 产量上浮 10%价格回落 5% 平年: (1.00, 1.00), 歉年: (0.85, 1.08)} # 产量下降 15%价格上浮 8% prob {丰年: 0.3, 平年: 0.5, 歉年: 0.2} total 0.0 for t in range(years): for name, (y_f, p_f) in scenarios.items(): profit_t 0.0 for i in range(len(land)): for j in range(len(crop)): profit_t ( X[t, i, j] * crop.loc[j, 亩产量] * y_f * crop.loc[j, 销售单价] * p_f - X[t, i, j] * crop.loc[j, 种植成本] ) total prob[name] * profit_t return -total penalty(X, land) # 负号让 PSO 从最大化转最小化expected_profit的目标值是对七个年份、三个情景的加权求和。penalty(X, land)是面积越界和轮作约束的罚函数实现我通常把越界面积平方乘一个大系数加回去而不是直接截断变量因为截断会破坏 PSO 搜索的连续性。3.2 粒子群参数怎么设置pso PSO(funcexpected_profit, n_dimyears * len(land) * len(crop), pop40, max_iter200, lb0, ubland[面积].max(), w0.8, c11.4, c21.4) pso.run() best_x pso.gbest_x这里的pop40是粒子数max_iter200是最大迭代轮数lb0和ubland[面积].max()给出变量上下界。w0.8是惯性权重控制粒子保持原有速度的程度c11.4是自我认知系数c21.4是社会认知系数两者相等会让粒子在局部搜索和全局收敛之间保持平衡。如果赛题时间紧张可以把max_iter降到 100但代价是结果不稳定。对比线性规划的单次确定解PSO 要跑 200 × 40 次目标函数评估Q2.py 运行十几分钟是正常的。与 Lingo 或 Matlab 里的ga工具箱相比sko.PSO不用写完整的约束矩阵只需要把硬约束塞进目标函数。代价是每次运行结果有随机性我建议运行三次取最优同时记录每次的期望利润浮动范围作为模型鲁棒性的参考指标。4. 销售量预测链路pearson.py、Linear_Regression.py 与 nonlinear_regression.py 的组合4.1 先检验相关性再决定用哪种模型问题三需要预测每种作物的预期销售量但不同作物的销售量随时间呈现出完全不同的形态有的近似直线上升有的先增后平有的周期性波动。pearson.py的作用是在做回归之前先回答一个问题这组年份和销售量之间到底存不存在线性相关如果皮尔逊系数绝对值低于 0.8 或 p 值大于 0.05直接套线性回归就是无效建模必须进入非线性分支。from scipy.stats import pearsonr def check_linearity(data, crop_name, year_col年份, sales_col预期销售量): sub data[data[作物名称] crop_name] r, p pearsonr(sub[year_col], sub[sales_col]) # r 在 [-1, 1]p 小于 0.05 说明相关性显著 return r, p for crop_name in data[作物名称].unique(): r, p check_linearity(data, crop_name) branch linear if abs(r) 0.8 and p 0.05 else nonlinear print(f{crop_name}: r{r:.3f}, p{p:.3f}, - {branch})赛题场景里样本量通常很小一个作物往往只有几年数据p 值的参考意义要打折。我的处理方式是同时看 r 和残差图如果 r 大但数据点明显弯成弧形仍然归入非线性分支如果 r 中等但 p 值极小说明线性趋势存在但噪声占比高应该线性与非线性各跑一版对比。4.2 线性回归用 statsmodels 而不是 sklearnLinear_Regression.py做的是普通最小二乘回归但有一个容易被忽略的细节sklearn 的LinearRegression不输出系数显著性、R² 和 F 统计量用来预测没问题用来解释模型就缺了诊断信息。比赛论文里需要写「模型显著性强」所以代码里更适合用 statsmodels。import statsmodels.api as sm def fit_ols(x, y): X sm.add_constant(x) # 加截距项 model sm.OLS(y, X).fit() return model.params, model.rsquared, model.pvaluesadd_constant是必要步骤缺失时模型会被强行过原点这对销售量数据是错误假设。model.rsquared就是拟合优度论文里报这个值model.pvalues用于判断年份变量是否显著通常小于 0.05 才会保留。4.3 非线性回归的曲线族选择nonlinear_regression.py对应三种最常用的曲线形态from scipy.optimize import curve_fit def quadratic(x, a, b, c): return a * x**2 b * x c def exp_growth(x, a, b, c): return a * np.exp(b * (x - x[0])) c def log_growth(x, a, b, c): return a * np.log(b * (x - x[0]) 1) cquadratic适合先升后降或持续加速的作物exp_growth适合处于快速扩张期的品种log_growth适合市场趋于饱和、增长逐步放缓的品种。curve_fit返回的popt是最优参数数组pcov是参数协方差矩阵对角线开方就是参数标准误。如果拟合后某个参数的标准误比参数本身还大说明曲线选型有问题应该换一个函数再拟合而不是直接采用结果。线性回归与非线性回归之间没有绝对优劣关键在前期判断。pearson.py 负责把「能不能用线性」这个问题量化后两个脚本负责在各自分支里产出预测值这条链路先检验再建模的顺序比拿到数据直接 fit 更经得起复现和答辩。5. 地块与大棚的差异化处理Q3_cluster.py 的 KMeans 与 Q3_ABCD_EF.py 的手动二分类5.1 为什么问题三要先聚类再预测问题三的数据里每种作物的预期销售量、销售单价、种植成本各有各的波动节奏把所有作物放进同一个回归模型会造成严重的异方差问题比如大宗粮食作物数量级是蔬菜的几十倍一起拟合时蔬菜的预测误差会被粮食作物的残差稀释掉。Q3_cluster.py的思路是先按销售行为把作物分成三类让同类作物拥有相近的参数尺度再用类别数据单独做回归。因为资源附带了「相关性信息_类别1.xlsx」「相关性信息_类别2.xlsx」「相关性信息_类别3.xlsx」三个输出文件说明 Q3_cluster.py 实际使用的就是三分类。聚类特征选取的是预期销售量、种植成本、销售单价三个维度的均值向量。from sklearn.cluster import KMeans from sklearn.preprocessing import StandardScaler def cluster_crops(df, feature_cols, n_clusters3): X df[feature_cols].values X StandardScaler().fit_transform(X) # 消除量纲差异 km KMeans(n_clustersn_clusters, random_state42, n_init10) df[类别] km.fit_predict(X) return dfStandardScaler这一步非常关键。预期销售量的方差远大于种植成本如果不做标准化聚类结果基本只被销售量主导类别 2 和类别 3 的区别会变成单纯的「量级差异」而不是「行为差异」。n_init10是让 KMeans 用 10 个不同的随机中心点各跑一次取最优结果避免初始化导致的局部最优。类别数 3 不是随便定的可以用肘部法则验证画类别数从 2 到 6 的 SSE 曲线选拐点。5.2 类别特征如何对应到回归聚类完成后输出文件「相关性信息_类别x.xlsx」里每一类包含一组作物。接下来的回归不是对所有作物混在一起做而是每个类别单独执行第 4 章的 pearson 回归流程。类别之间的差异往往体现在曲线形态上类别 1 可能是量价齐升的经济作物类别 2 是价格平稳的粮食作物类别 3 是价格随季节波动大的果蔬类。这样分类后回归得到的 R² 会比整体拟合高出不少论文里说明「先聚类再分类预测」也就是在讲这一步。5.3 露天地块与大棚的拆分Q3_ABCD_EF.py是一个手动二分类脚本逻辑比 KMeans 更直接文件名里的 A/B/C/D 对应赛题中的四类露天地块E/F 对应大棚地块。大棚作物的生长环境可控销售量与温度的耦合程度低于露天作物混合建模会引入难以解释的干扰项。常见做法是先用地块类型字段把数据拆成「露天」和「大棚」两份再分别走预测流程。def split_land_category(land): outdoor [平旱地, 梯田, 山坡地, 水浇地] land[地块大类] land[地块类型].apply( lambda t: 露天 if t in outdoor else 大棚) return land[land[地块大类] 露天], land[land[地块大类] 大棚]这里apply配合 lambda 做映射比逐行for循环快一个数量级。拆分之后两份数据各自做一次聚类和回归输出对应为「result3(地块-大棚).xlsx」和「result3(聚类).xlsx」。两个结果文件不是替代关系前者回答地块层面的分配约束后者回答作物层面的销售量预测最终方案要把两者合并先用聚类结果算出每类作物的预测销售量再用地块-大棚分组把销售量落实到具体地块面积上。6. 环境复现、运行顺序与结果交叉验证6.1 依赖库安装与版本坑资源标注要求 Python 3.11核心依赖是 pandas、numpy、sko、matplotlib。其中sko是 scikit-opt 的缩写PyPI 包名是scikit-opt直接pip install sko会装错包。pip install pandas numpy scipy scikit-learn scikit-opt matplotlib openpyxlopenpyxl是 pandas 读写 .xlsx 的底层引擎不装会在pd.read_excel时报 ImportError。scikit-opt 目前对 numpy 2.x 的兼容性有历史问题建议先安装 numpy1.26.4 再装 scikit-opt避免出现module numpy has no attribute float这类报错。statsmodels 只在回归脚本里用到如果只跑 Q3 的线性回归分支可以用 sklearn 替代但建议还是装上。6.2 脚本执行顺序顺序脚本输入输出1Q1_1.py地块信息.xlsx、作物信息.xlsxresult1_1.xlsx2Q1_2.py地块信息.xlsx、作物信息.xlsxresult1_2.xlsx3Q2.py地块信息.xlsx、作物信息.xlsx、汇总补充信息.xlsxresult2.xlsx4Q3_cluster.py相关性信息.xlsx相关性信息_类别1/2/3.xlsx5pearson.py相关性信息_类别1/2/3.xlsx相关系数与分支标记6Linear_Regression.py nonlinear_regression.py类别数据预期销售量预测值7Q3_ABCD_EF.py地块信息.xlsx、作物信息.xlsx、预测值result3(地块-大棚).xlsx、result3(聚类).xlsx6.3 一个值得养成的验证技巧def verify_area_constraint(solution, land, eps1e-6): # solution 是 predict 或 优化结果行是地块列是作物 total_planted solution.sum(axis1) # 每块地的总种植面积 exceed total_planted - land[面积] return exceed[exceed eps]这套资源的问题在于脚本多、中间文件多最容易出的错不是模型错误而是某个中间文件被上一步覆盖导致后续脚本读到错位数据。我的习惯是跑完 Q3_cluster.py 和 Q3_ABCD_EF.py 之后第一时间用上面的函数检查「result3(聚类).xlsx」里每块地的种植面积总和是否仍然满足约束。如果超过问题大概率出在聚类后的作物归属发生了错位如果没超过但方案利润异常低再回去查回归脚本里是不是把某个类别的预测销售量串到了另一个类别。复现这套代码时把握好「输入文件路径正确」「输出文件不互相覆盖」「约束校验通过」这三条整套流程基本不会出问题。本文还有配套的精品资源点击获取
数学建模C题Python实现：从线性规划到粒子群与聚类预测
简介针对2024年全国大学生数学建模竞赛C题这份资源围绕农作物种植策略提供了一整套可落地的分析方案涵盖建模思路、程序源码和优秀论文适合参赛学生、赛事指导老师以及想入门数学建模的读者使用。代码按赛题问题拆分前两问均有独立求解脚本第三问综合运用线性回归、非线性回归、Pearson相关性检验和KMeans聚类算法对预期销售量、种植成本和销售单价进行统计建模并对地块—大棚做分类处理还附带整理好的地块信息、作物信息、相关性信息及2023年汇总补充数据基本还原了从数据处理到结果输出的完整链路。资源共28个文件包括12个Excel数据文件、8个Python程序脚本、2个PDF赛题与论文、2个Markdown说明文档及4个代码备份文件压缩包整体仅2.21MB轻量易用方便按需查阅。当前已有133人学习适合需要复现国赛获奖思路、参考建模代码或快速了解农作物种植优化问题的读者。1. 从赛题到代码一份完整 C 题方案的解构2024 年国赛 C 题把「农作物种植策略」做成了一个相当实际的优化问题村里有平旱地、梯田、山坡地、水浇地还有大棚十几种作物要在 2024—2030 年间逐年排产产量、售价、成本还会波动。这份资源里真正值钱的不是那篇论文而是六个按问题顺序编号的 Python 求解文件从 Q1_1.py、Q1_2.py 的线性规划起步到 Q2.py 的随机优化再到问题三的 pearson.py、线性/非线性回归脚本和 KMeans 聚类完整覆盖了确定性优化、不确定性建模、相关性检验和聚类四类数学建模模型。它适合两类人正在备战数学建模竞赛、想找一套可复现代码路径的选手以及用 Python 处理生产调度、想把优化工具链串起来落地的开发者。下文按赛题推进顺序逐层拆这套代码。2. 地块与作物的约束系统Q1_1.py 和 Q1_2.py 的线性规划实现2.1 从赛题数据到决策变量问题一的核心假设是未来产量、售价、成本相对稳定因此方案可以在同一套约束下逐年复用。数据基础是资源包里的「地块信息.xlsx」和「作物信息.xlsx」前者记录每块地的类型与面积后者记录每种作物的亩产量、种植成本、销售单价。建模的第一步是把自然语言描述转成线性规划的标准形式。决策变量定义为# x[i * n_crop j] 表示第 i 块地上种植第 j 种作物的面积 # 索引展开方式地块优先同一地块内的作物连续排列 n_land len(land) n_crop len(crop) x [0] * (n_land * n_crop)变量采用「地块优先 块内作物连续」的扁平排列是因为scipy.optimize.linprog只接受一维决策变量后续构造约束矩阵时需要按i * n_crop j这个索引规则定位。如果换成二维变量矩阵构造会直观一些但还原成线性规划标准形式时需要多一层 reshape反而容易出错。2.2 目标函数与约束矩阵的构造目标函数是总利润最大化。单块地上某作物的利润按「面积 × 亩产量 × 销售单价 − 面积 × 种植成本」计算因此目标系数向量里的每一项都是单位面积净利润。linprog只做最小化所以把利润向量整体取负。import pandas as pd from scipy.optimize import linprog def build_q1_model(land, crop): n_land, n_crop len(land), len(crop) c [] for _, ld in land.iterrows(): for _, cr in crop.iterrows(): profit cr[亩产量] * cr[销售单价] - cr[种植成本] c.append(-profit) # 最小化负利润 最大化利润 A_ub, b_ub [], [] # 约束每块地上所有作物的种植面积之和 该地块面积 for i in range(n_land): row [0] * (n_land * n_crop) for j in range(n_crop): row[i * n_crop j] 1 A_ub.append(row) b_ub.append(land.loc[i, 面积]) # 变量下界为 0上界为该地块面积保证单个变量不会超过地块总面积 bounds [(0, land.loc[i, 面积]) for i in range(n_land) for _ in range(n_crop)] res linprog(c, A_ubA_ub, b_ubb_ub, boundsbounds, methodhighs) return res这里的methodhighs是 SciPy 1.7 之后默认启用的 HiGHS 求解器适合中等规模的线性规划如果直接用旧代码里的methodsimplex在数千个变量的场景下迭代步数和数值误差都会明显变大。约束矩阵 A_ub 的每一行只在一块地对应的 n_crop 个位置上填 1本质是块内加和约束这在实际求解中比单独对每块地每种作物分别设限更接近赛题的「可以少种但不可以多种」逻辑。2.3 Q1_2.py 在 Q1_1 之上叠加的时间维约束Q1_1.py 求解的是稳定假设下的一期最优布局而 Q1_2.py 解决的是问题一的第二问方案要逐年滚动到 2030 年同一块地上不能连续多年种植同一种作物。此时决策变量需要扩成三维x[t][i][j]表示第 t 年第 i 块地种第 j 种作物的面积并额外追加重茬约束。T 7 # 2024 到 2030 def add_rotation_constraint(A_ub, b_ub, n_land, n_crop, T, forbid_pairs): for t in range(1, T): # 从第二年开始检查上一年的种植情况 for i in range(n_land): for j in forbid_pairs.get(i, []): row [0] * (T * n_land * n_crop) # 上一年的面积 这一年的面积 上一年面积的阈值 row[(t-1) * n_land * n_crop i * n_crop j] 1 row[t * n_land * n_crop i * n_crop j] 0.2 A_ub.append(row) b_ub.append(0) return A_ub, b_ub这段约束的含义是如果第 t−1 年某地块种了易重茬障碍的作物第 t 年的同作物面积不能超过上一年的 0.2 倍相当于允许少量留茬但不允许连续性大规模种植。forbid_pairs需要根据作物信息里的「重茬敏感」标记手动配置不同作物对连作的容忍度差异很大这属于赛题数据里不会直接给出、但实际求解时必须补全的工程细节。Q1_1 和 Q1_2 的关系可以简单理解成前者是「一年定终身」后者是「七年滚动排产」。从代码量看Q1_2 只多了时间索引和轮作约束但求解规模扩大了近 7 倍这也是我建议用highs而不是老求解器的直接原因。3. 不确定性建模Q2.py 与改进粒子群优化算法IPSO3.1 为什么线性规划在问题二失效问题二把产量、售价、成本同时变成随机变量。如果继续用线性规划就得把每个随机参数替换成期望值但期望值解在真实波动下往往违背面积约束或产生明显亏损。原因是线性规划天然假设参数确定一旦目标函数里的系数带随机性最优解对参数扰动非常敏感尤其是产量和价格呈负相关时丰年产量高但价格低直接用均值建模会系统性高估收益。资源包里的论文标题给出了作者团队的解法基于粒子群优化算法的华北农作物种植策略模型。代码里对应的是 Q2.py它引入三个情景——丰年、平年、歉年每个情景有独立的产量折减系数和价格波动系数目标函数改为多年期望利润。import numpy as np from sko.PSO import PSO def expected_profit(x, land, crop, years7): X x.reshape(years, len(land), len(crop)) scenarios {丰年: (1.10, 0.95), # 产量上浮 10%价格回落 5% 平年: (1.00, 1.00), 歉年: (0.85, 1.08)} # 产量下降 15%价格上浮 8% prob {丰年: 0.3, 平年: 0.5, 歉年: 0.2} total 0.0 for t in range(years): for name, (y_f, p_f) in scenarios.items(): profit_t 0.0 for i in range(len(land)): for j in range(len(crop)): profit_t ( X[t, i, j] * crop.loc[j, 亩产量] * y_f * crop.loc[j, 销售单价] * p_f - X[t, i, j] * crop.loc[j, 种植成本] ) total prob[name] * profit_t return -total penalty(X, land) # 负号让 PSO 从最大化转最小化expected_profit的目标值是对七个年份、三个情景的加权求和。penalty(X, land)是面积越界和轮作约束的罚函数实现我通常把越界面积平方乘一个大系数加回去而不是直接截断变量因为截断会破坏 PSO 搜索的连续性。3.2 粒子群参数怎么设置pso PSO(funcexpected_profit, n_dimyears * len(land) * len(crop), pop40, max_iter200, lb0, ubland[面积].max(), w0.8, c11.4, c21.4) pso.run() best_x pso.gbest_x这里的pop40是粒子数max_iter200是最大迭代轮数lb0和ubland[面积].max()给出变量上下界。w0.8是惯性权重控制粒子保持原有速度的程度c11.4是自我认知系数c21.4是社会认知系数两者相等会让粒子在局部搜索和全局收敛之间保持平衡。如果赛题时间紧张可以把max_iter降到 100但代价是结果不稳定。对比线性规划的单次确定解PSO 要跑 200 × 40 次目标函数评估Q2.py 运行十几分钟是正常的。与 Lingo 或 Matlab 里的ga工具箱相比sko.PSO不用写完整的约束矩阵只需要把硬约束塞进目标函数。代价是每次运行结果有随机性我建议运行三次取最优同时记录每次的期望利润浮动范围作为模型鲁棒性的参考指标。4. 销售量预测链路pearson.py、Linear_Regression.py 与 nonlinear_regression.py 的组合4.1 先检验相关性再决定用哪种模型问题三需要预测每种作物的预期销售量但不同作物的销售量随时间呈现出完全不同的形态有的近似直线上升有的先增后平有的周期性波动。pearson.py的作用是在做回归之前先回答一个问题这组年份和销售量之间到底存不存在线性相关如果皮尔逊系数绝对值低于 0.8 或 p 值大于 0.05直接套线性回归就是无效建模必须进入非线性分支。from scipy.stats import pearsonr def check_linearity(data, crop_name, year_col年份, sales_col预期销售量): sub data[data[作物名称] crop_name] r, p pearsonr(sub[year_col], sub[sales_col]) # r 在 [-1, 1]p 小于 0.05 说明相关性显著 return r, p for crop_name in data[作物名称].unique(): r, p check_linearity(data, crop_name) branch linear if abs(r) 0.8 and p 0.05 else nonlinear print(f{crop_name}: r{r:.3f}, p{p:.3f}, - {branch})赛题场景里样本量通常很小一个作物往往只有几年数据p 值的参考意义要打折。我的处理方式是同时看 r 和残差图如果 r 大但数据点明显弯成弧形仍然归入非线性分支如果 r 中等但 p 值极小说明线性趋势存在但噪声占比高应该线性与非线性各跑一版对比。4.2 线性回归用 statsmodels 而不是 sklearnLinear_Regression.py做的是普通最小二乘回归但有一个容易被忽略的细节sklearn 的LinearRegression不输出系数显著性、R² 和 F 统计量用来预测没问题用来解释模型就缺了诊断信息。比赛论文里需要写「模型显著性强」所以代码里更适合用 statsmodels。import statsmodels.api as sm def fit_ols(x, y): X sm.add_constant(x) # 加截距项 model sm.OLS(y, X).fit() return model.params, model.rsquared, model.pvaluesadd_constant是必要步骤缺失时模型会被强行过原点这对销售量数据是错误假设。model.rsquared就是拟合优度论文里报这个值model.pvalues用于判断年份变量是否显著通常小于 0.05 才会保留。4.3 非线性回归的曲线族选择nonlinear_regression.py对应三种最常用的曲线形态from scipy.optimize import curve_fit def quadratic(x, a, b, c): return a * x**2 b * x c def exp_growth(x, a, b, c): return a * np.exp(b * (x - x[0])) c def log_growth(x, a, b, c): return a * np.log(b * (x - x[0]) 1) cquadratic适合先升后降或持续加速的作物exp_growth适合处于快速扩张期的品种log_growth适合市场趋于饱和、增长逐步放缓的品种。curve_fit返回的popt是最优参数数组pcov是参数协方差矩阵对角线开方就是参数标准误。如果拟合后某个参数的标准误比参数本身还大说明曲线选型有问题应该换一个函数再拟合而不是直接采用结果。线性回归与非线性回归之间没有绝对优劣关键在前期判断。pearson.py 负责把「能不能用线性」这个问题量化后两个脚本负责在各自分支里产出预测值这条链路先检验再建模的顺序比拿到数据直接 fit 更经得起复现和答辩。5. 地块与大棚的差异化处理Q3_cluster.py 的 KMeans 与 Q3_ABCD_EF.py 的手动二分类5.1 为什么问题三要先聚类再预测问题三的数据里每种作物的预期销售量、销售单价、种植成本各有各的波动节奏把所有作物放进同一个回归模型会造成严重的异方差问题比如大宗粮食作物数量级是蔬菜的几十倍一起拟合时蔬菜的预测误差会被粮食作物的残差稀释掉。Q3_cluster.py的思路是先按销售行为把作物分成三类让同类作物拥有相近的参数尺度再用类别数据单独做回归。因为资源附带了「相关性信息_类别1.xlsx」「相关性信息_类别2.xlsx」「相关性信息_类别3.xlsx」三个输出文件说明 Q3_cluster.py 实际使用的就是三分类。聚类特征选取的是预期销售量、种植成本、销售单价三个维度的均值向量。from sklearn.cluster import KMeans from sklearn.preprocessing import StandardScaler def cluster_crops(df, feature_cols, n_clusters3): X df[feature_cols].values X StandardScaler().fit_transform(X) # 消除量纲差异 km KMeans(n_clustersn_clusters, random_state42, n_init10) df[类别] km.fit_predict(X) return dfStandardScaler这一步非常关键。预期销售量的方差远大于种植成本如果不做标准化聚类结果基本只被销售量主导类别 2 和类别 3 的区别会变成单纯的「量级差异」而不是「行为差异」。n_init10是让 KMeans 用 10 个不同的随机中心点各跑一次取最优结果避免初始化导致的局部最优。类别数 3 不是随便定的可以用肘部法则验证画类别数从 2 到 6 的 SSE 曲线选拐点。5.2 类别特征如何对应到回归聚类完成后输出文件「相关性信息_类别x.xlsx」里每一类包含一组作物。接下来的回归不是对所有作物混在一起做而是每个类别单独执行第 4 章的 pearson 回归流程。类别之间的差异往往体现在曲线形态上类别 1 可能是量价齐升的经济作物类别 2 是价格平稳的粮食作物类别 3 是价格随季节波动大的果蔬类。这样分类后回归得到的 R² 会比整体拟合高出不少论文里说明「先聚类再分类预测」也就是在讲这一步。5.3 露天地块与大棚的拆分Q3_ABCD_EF.py是一个手动二分类脚本逻辑比 KMeans 更直接文件名里的 A/B/C/D 对应赛题中的四类露天地块E/F 对应大棚地块。大棚作物的生长环境可控销售量与温度的耦合程度低于露天作物混合建模会引入难以解释的干扰项。常见做法是先用地块类型字段把数据拆成「露天」和「大棚」两份再分别走预测流程。def split_land_category(land): outdoor [平旱地, 梯田, 山坡地, 水浇地] land[地块大类] land[地块类型].apply( lambda t: 露天 if t in outdoor else 大棚) return land[land[地块大类] 露天], land[land[地块大类] 大棚]这里apply配合 lambda 做映射比逐行for循环快一个数量级。拆分之后两份数据各自做一次聚类和回归输出对应为「result3(地块-大棚).xlsx」和「result3(聚类).xlsx」。两个结果文件不是替代关系前者回答地块层面的分配约束后者回答作物层面的销售量预测最终方案要把两者合并先用聚类结果算出每类作物的预测销售量再用地块-大棚分组把销售量落实到具体地块面积上。6. 环境复现、运行顺序与结果交叉验证6.1 依赖库安装与版本坑资源标注要求 Python 3.11核心依赖是 pandas、numpy、sko、matplotlib。其中sko是 scikit-opt 的缩写PyPI 包名是scikit-opt直接pip install sko会装错包。pip install pandas numpy scipy scikit-learn scikit-opt matplotlib openpyxlopenpyxl是 pandas 读写 .xlsx 的底层引擎不装会在pd.read_excel时报 ImportError。scikit-opt 目前对 numpy 2.x 的兼容性有历史问题建议先安装 numpy1.26.4 再装 scikit-opt避免出现module numpy has no attribute float这类报错。statsmodels 只在回归脚本里用到如果只跑 Q3 的线性回归分支可以用 sklearn 替代但建议还是装上。6.2 脚本执行顺序顺序脚本输入输出1Q1_1.py地块信息.xlsx、作物信息.xlsxresult1_1.xlsx2Q1_2.py地块信息.xlsx、作物信息.xlsxresult1_2.xlsx3Q2.py地块信息.xlsx、作物信息.xlsx、汇总补充信息.xlsxresult2.xlsx4Q3_cluster.py相关性信息.xlsx相关性信息_类别1/2/3.xlsx5pearson.py相关性信息_类别1/2/3.xlsx相关系数与分支标记6Linear_Regression.py nonlinear_regression.py类别数据预期销售量预测值7Q3_ABCD_EF.py地块信息.xlsx、作物信息.xlsx、预测值result3(地块-大棚).xlsx、result3(聚类).xlsx6.3 一个值得养成的验证技巧def verify_area_constraint(solution, land, eps1e-6): # solution 是 predict 或 优化结果行是地块列是作物 total_planted solution.sum(axis1) # 每块地的总种植面积 exceed total_planted - land[面积] return exceed[exceed eps]这套资源的问题在于脚本多、中间文件多最容易出的错不是模型错误而是某个中间文件被上一步覆盖导致后续脚本读到错位数据。我的习惯是跑完 Q3_cluster.py 和 Q3_ABCD_EF.py 之后第一时间用上面的函数检查「result3(聚类).xlsx」里每块地的种植面积总和是否仍然满足约束。如果超过问题大概率出在聚类后的作物归属发生了错位如果没超过但方案利润异常低再回去查回归脚本里是不是把某个类别的预测销售量串到了另一个类别。复现这套代码时把握好「输入文件路径正确」「输出文件不互相覆盖」「约束校验通过」这三条整套流程基本不会出问题。本文还有配套的精品资源点击获取
