ArcGIS Pro随机森林预测实战：从数据准备到模型部署全流程解析

1次阅读

没有评论

共计 3762 个字符，预计需要花费 10 分钟才能阅读完成。

随机森林作为一种集成学习方法，在 ArcGIS Pro 中特别适合处理具有空间特性的预测任务。以下是几个典型应用场景：

土地利用 / 土地覆盖分类：利用多光谱影像和其他辅助数据（如 DEM）进行高精度分类
自然灾害风险评估：基于历史灾害数据预测山火、滑坡等事件的发生概率
城市扩张模拟：结合社会经济因子预测未来城市发展热点区域
环境指标预测：如空气质量、地下水污染等空间连续变量的插值

这些场景的共同特点是：数据具有空间自相关性、包含多源异构特征、且需要处理非线性关系。

特征工程复杂：空间数据通常包含坐标、邻域关系等特殊特征
样本不平衡：某些类别（如灾害发生点）在样本中占比极低
计算效率低下：传统方法难以处理大规模栅格数据

使用 ArcPy 进行空间特征提取（如缓冲区统计、地形指标计算）
采用分层抽样解决类别不平衡问题
通过并行计算和内存优化提升处理效率

import arcpy
from sklearn.ensemble import RandomForestClassifier

# 设置工作空间
arcpy.env.workspace = "C:/Data/Project.gdb"
arcpy.env.overwriteOutput = True

# 坐标系转换（示例：WGS84 转 UTM）input_fc = "raw_data"
output_fc = "projected_data"
out_coordinate_system = arcpy.SpatialReference(32651)  # WGS84 UTM Zone 51N
arcpy.Project_management(input_fc, output_fc, out_coordinate_system)

# 处理空值（以 NDVI 字段为例）with arcpy.da.UpdateCursor(output_fc, ["NDVI"]) as cursor:
    for row in cursor:
        if row[0] is None:
            row[0] = 0  # 用 0 填充缺失值
            cursor.updateRow(row)

# 计算邻域统计特征（以 500 米缓冲区为例）arcpy.Statistics_analysis("projected_data", "buffer_stats", 
                         [["Elevation", "MEAN"], ["Slope", "MAX"]], 
                         "NEIGHBORHOOD", "Circle 500 METERS")

# 创建空间权重矩阵（用于检测自相关性）arcpy.GenerateSpatialWeightsMatrix_stats("projected_data", "ID",
                                        "spatial_weights", "K_NEAREST_NEIGHBORS",
                                        "EUCLIDEAN", 12)  # 选择 12 个最近邻

import pandas as pd
from sklearn.model_selection import train_test_split

# 读取属性表
arr = arcpy.da.TableToNumPyArray("final_features", ["*"])
df = pd.DataFrame(arr)

# 划分训练测试集（考虑空间分块）X = df.drop(["CLASS", "Shape"], axis=1)
y = df["CLASS"]
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=42)

# 初始化随机森林（启用 OOB 评估）rf = RandomForestClassifier(n_estimators=500,
                           max_depth=15,
                           oob_score=True,
                           n_jobs=-1,  # 使用所有 CPU 核心
                           class_weight="balanced")  # 处理类别不平衡

# 训练模型
rf.fit(X_train, y_train)

# 评估指标
print(f"OOB Score: {rf.oob_score_:.3f}")
print(f"Test Accuracy: {rf.score(X_test, y_test):.3f}")

# 特征重要性可视化
import matplotlib.pyplot as plt
plt.barh(X.columns, rf.feature_importances_)
plt.xlabel("Feature Importance")
plt.show()

n_jobs=-1：使用所有可用 CPU 核心
max_samples=0.8：每棵树使用 80% 样本加快训练
bootstrap=True：启用自助采样减少过拟合

# 分块处理大栅格数据
arcpy.env.compression = "LZW"
arcpy.env.pyramid = "PYRAMIDS -1"
arcpy.env.tileSize = "128 128"  # 设置分块大小

# 增量训练（适用于新增数据）from sklearn.ensemble import RandomForestClassifier
rf = RandomForestClassifier(warm_start=True)  # 允许增量训练

# 首次训练
rf.fit(X_train_part1, y_train_part1)

# 后续追加数据
rf.n_estimators += 100
rf.fit(X_train_part2, y_train_part2)

import joblib

# 保存模型
joblib.dump(rf, "rf_model.pkl")

# 在 ArcPy 中加载使用
loaded_model = joblib.load("rf_model.pkl")
prediction = loaded_model.predict_proba(new_data)[:, 1]

# 将预测结果写回要素类
with arcpy.da.UpdateCursor("prediction_points", ["PROBABILITY"]) as cursor:
    for i, row in enumerate(cursor):
        row[0] = prediction[i]
        cursor.updateRow(row)

解决方法：
设置 class_weight="balanced" 参数
采用分层抽样（train_test_split的 stratify 参数）
使用 SMOTE 过采样技术

Moran’s I 检验：

from pysal.explore import esda
import numpy as np

# 计算空间自相关
w = ... # 空间权重矩阵
moran = esda.Moran(y_train, w)
print(f"Moran's I: {moran.I}, p-value: {moran.p_sim}")

# 若存在显著自相关（p<0.05）：# 1. 加入空间滞后变量作为特征
# 2. 采用空间交叉验证

距离相关特征必须使用投影坐标系
建议选择 UTM 或 Albers 等面积投影
通过 arcpy.Describe().spatialReference 检查当前坐标系

指标	随机森林	逻辑回归	克里金插值
处理非线性关系	优秀	需手动添加交互项	中等
特征重要性	内置评估	需要额外计算	不适用
计算效率	并行化优势	高效	大区域较慢

创建 Python 脚本工具：

import arcpy

class Toolbox(object):
    def __init__(self):
        self.label = "Random Forest Predictor"
        self.alias = ""
        self.tools = [PredictTool]

class PredictTool(object):
    def __init__(self):
        self.label = "Predict using RF"
        self.description = ""
        self.canRunInBackground = True

    def execute(self, parameters, messages):
        input_fc = parameters[0].valueAsText
        model_path = parameters[1].valueAsText
        output_fc = parameters[2].valueAsText

        # 调用前面实现的预测逻辑
        arcpy.AddMessage("Prediction completed!")
        return