# 城市创新驱动因素可解释性机器学习分析指南

本文档包含了复现本研究中 XGBoost 与 SHAP 分析流程的详细步骤和运行指令。所有的分析代码均保存在 `数据和代码` 文件夹下。

---

## 🛠 运行环境准备

在运行代码之前，请确保你的 Python 环境中已经安装了以下依赖库。你可以在终端中运行以下指令进行安装：

```bash
pip install pandas numpy scikit-learn xgboost lightgbm catboost shap matplotlib
```

---

## 🚀 脚本执行顺序与说明

请**严格按照以下顺序**在终端中运行脚本，以确保数据流的正确传递。

### 1. 面板数据构建与合并
- **脚本文件**：`build_panel.py`
- **功能**：从原始的 Excel/DTA/CSV 格式源文件中，提取 102 个目标城市的数据，基于“协作启动年份”进行面板数据截断，并使用**熵权法 (EWM)**计算得到因变量 `Innovation_Index`（综合创新指数）。
- **运行指令**：
  ```bash
  python build_panel.py
  ```
- **输出**：生成基础的面板数据集 `面板数据最终版.csv`。

---

### 2. 数据缺失值高级插补
- **脚本文件**：`data_imputation.py`
- **功能**：针对部分西部城市存在的指标缺失（如外贸依存度、互联网普及率等），采用“面板内双向线性插值”与“基于 XGBoost 的多重迭代插补（IterativeImputer）”相结合的两步法进行处理。
- **运行指令**：
  ```bash
  python data_imputation.py
  ```
- **输出**：生成零缺失的高质量数据集 `面板数据_多重插补后.csv`，作为后续所有机器学习模型的统一输入。

---

### 3. 模型超参数优化 (可选)
- **脚本文件**：`hyperparameter_tuning.py`
- **功能**：使用 5 折交叉验证（5-fold CV）和 RandomizedSearchCV 对 XGBoost 的超参数进行网格搜索。旨在找到能够兼顾预测精度与特征利用率的最佳参数组合（特别是调整了列采样率以暴露边缘特征）。
- **运行指令**：
  ```bash
  python hyperparameter_tuning.py
  ```
- **输出**：在终端打印出最佳超参数字典及 R² 得分。

---

### 4. 多算法性能对比
- **脚本文件**：`model_comparison.py`
- **功能**：对齐论文 Table 3 的研究设计，对比 XGBoost, LightGBM, CatBoost, Random Forest, GBDT, SVR 以及传统 OLS 的模型表现。证明树集成模型在处理非线性创新驱动关系时的绝对优势。
- **运行指令**：
  ```bash
  python model_comparison.py
  ```
- **输出**：在 `模型结果与图片` 文件夹下生成对比表格 `Model_Comparison_Table_3.csv`。

---

### 5. SHAP 与 PDP 可解释性分析 (核心)
- **脚本文件**：`xgboost_shap_analysis.py`
- **功能**：本研究的核心分析脚本。它将训练好的 XGBoost 模型输入到 SHAP 解释器中，进行全局特征重要性归因，并计算部分依赖图 (PDP) 以寻找阈值效应。
- **运行指令**：
  ```bash
  python xgboost_shap_analysis.py
  ```
- **输出**：在 `模型结果与图片` 文件夹下生成以下论文级图表：
  1. `Feature_Importance_Table_Innovation_Index.csv`：各特征的贡献百分比与正负向影响表（对齐论文 Table 5）。
  2. `SHAP_Gradient_Bar_Importance_Innovation_Index.png`：带渐变色的全特征重要性排版图（对齐论文 Fig 4-7）。
  3. `PDP_Custom_Style_All_Innovation_Index.png`：包含 28 个子图、置信区间阴影、基准线的高级 PDP 阈值效应矩阵图。
  4. `SHAP_Summary_Plot_Innovation_Index.png`：展示具体样本分布的 SHAP 蜂群图。

---

## 💡 如何替换因变量？
目前的脚本默认以 `Innovation_Index`（综合创新指数）作为研究目标。如果你需要分析其他维度的创新表现（例如：`pat_inv` 发明专利数量，或者 `newHT` 新增高新企业数），只需执行以下操作：
1. 打开 `xgboost_shap_analysis.py` 或其他脚本。
2. 找到代码顶部的 `target = 'Innovation_Index'`。
3. 将其修改为对应的列名（如 `target = 'pat_inv'`）。
4. 重新运行对应的 Python 脚本，新的图表和表格会自动生成并加上该变量的后缀名，不会覆盖旧文件。