Python scikit-learn 1.4 实战：蒙特西尼奥森林火灾数据 K-Means 聚类与 3D 可视化

当数据分析遇上环境科学，Python 的 scikit-learn 库为我们打开了一扇探索自然现象的新窗口。蒙特西尼奥自然公园的火灾数据集就像一本等待破译的自然密码本，而 K-Means 算法和 3D 可视化技术则是我们解读这些密码的利器。本文将带你深入这个数据科学项目，从数据准备到最终可视化，一步步揭开隐藏在气象数据背后的森林火灾模式。

1. 环境准备与数据加载

在开始我们的分析之前，确保你的 Python 环境已经安装了以下必要库：

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.cluster import KMeans
from mpl_toolkits.mplot3d import Axes3D
from sklearn.preprocessing import StandardScaler

数据集可以从阿里云天池平台获取，包含了葡萄牙蒙特西尼奥自然公园 517 起火灾事件的详细记录。让我们先加载并初步探索数据：

# 加载数据集
data = pd.read_csv('montesinho_fire_dataset.csv')

# 查看数据概览
print(data.info())
print(data.describe())

数据集包含以下关键特征：

• FFMC ：细小可燃物湿度码
• DMC ：深层可燃物湿度码
• DC ：干旱码
• ISI ：初始蔓延指数
• Temp ：温度(℃)
• RH ：相对湿度(%)
• Wind ：风速(km/h)
• Rain ：降雨量(mm)
• Area ：过火面积(公顷)

提示：在实际项目中，建议将数据分为训练集和测试集，但聚类分析通常使用全部数据进行模式发现。

2. 数据预处理与特征工程

原始数据需要经过适当处理才能发挥最大价值。我们重点关注与火灾风险直接相关的五个关键指标：

# 选择关键特征
features = ['FFMC', 'DMC', 'DC', 'ISI', 'Temp']
X = data[features]

# 数据标准化
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

# 检查缺失值
print("缺失值统计：\n", X.isnull().sum())

标准化处理对 K-Means 算法尤为重要，因为该算法对特征的尺度敏感。我们使用 StandardScaler 将各特征缩放到均值为0，方差为1的标准正态分布。

特征相关性分析表 ：

特征组合	相关系数	物理意义
FFMC-DMC	0.62	地表与深层可燃物湿度关系
DMC-DC	0.78	深层可燃物与干旱程度关联
ISI-Temp	0.55	火势蔓延与温度相关性

3. K-Means 聚类实战

确定最佳聚类数是 K-Means 算法的关键挑战。我们使用肘部法则和轮廓系数相结合的方法：

from sklearn.metrics import silhouette_score

inertia = []
silhouette_scores = []
K_range = range(2, 10)

for k in K_range:
    kmeans = KMeans(n_clusters=k, random_state=42, n_init=10)
    kmeans.fit(X_scaled)
    inertia.append(kmeans.inertia_)
    silhouette_scores.append(silhouette_score(X_scaled, kmeans.labels_))

# 绘制肘部法则图
plt.figure(figsize=(12, 5))
plt.subplot(1, 2, 1)
plt.plot(K_range, inertia, 'bo-')
plt.xlabel('聚类数 k')
plt.ylabel('惯性值')
plt.title('肘部法则')

plt.subplot(1, 2, 2)
plt.plot(K_range, silhouette_scores, 'ro-')
plt.xlabel('聚类数 k')
plt.ylabel('轮廓系数')
plt.title('轮廓系数分析')
plt.show()

基于分析结果，我们选择 k=3 进行最终聚类：

# 最终聚类模型
optimal_k = 3
kmeans = KMeans(n_clusters=optimal_k, random_state=42, n_init=10)
clusters = kmeans.fit_predict(X_scaled)

# 将聚类结果添加到原始数据
data['Cluster'] = clusters

聚类中心特征对比 ：

特征	集群0	集群1	集群2
FFMC	86.3	91.2	82.1
DMC	35.6	108.4	15.2
DC	220.7	540.3	110.5
ISI	7.2	14.5	4.3
Temp	18.6	25.3	12.7

4. 3D 可视化实现

现在是时候让数据"立体"起来了。我们将使用 Matplotlib 的 3D 功能展示聚类结果：

# 创建3D图形
fig = plt.figure(figsize=(12, 10))
ax = fig.add_subplot(111, projection='3d')

# 为每个集群设置不同颜色
colors = ['r', 'g', 'b']
for i in range(optimal_k):
    cluster_data = data[data['Cluster'] == i]
    ax.scatter(cluster_data['FFMC'], 
               cluster_data['DMC'], 
               cluster_data['Temp'],
               c=colors[i],
               label=f'集群 {i}',
               s=50,
               alpha=0.6)

# 标记聚类中心
centers = scaler.inverse_transform(kmeans.cluster_centers_)
ax.scatter(centers[:, 0], centers[:, 1], centers[:, 4],
           c='black', marker='X', s=200, alpha=1, label='聚类中心')

# 设置坐标轴标签
ax.set_xlabel('FFMC (细小可燃物湿度)')
ax.set_ylabel('DMC (深层可燃物湿度)')
ax.set_zlabel('Temperature (℃)')
ax.set_title('森林火灾数据三维聚类可视化')
plt.legend()
plt.tight_layout()

# 添加交互式旋转功能（在Jupyter中启用）
%matplotlib notebook  
plt.show()

可视化优化技巧 ：

• 调整 alpha 参数控制点透明度，避免重叠点遮挡
• 使用 s 参数调整点的大小，突出重要数据
• 添加 %matplotlib notebook 魔法命令实现交互式旋转

5. 聚类结果分析与物理解释

三个集群展现出明显的环境特征差异：

1. 高风险集群（集群1） ：

   • 特征：高温度(25.3℃)、高DMC(108.4)、高DC(540.3)
   • 物理解释：代表极端干燥和高温条件，火灾发生风险最高
   • 实际意义：需要最高级别的防火警戒

2. 中风险集群（集群0） ：

   • 特征：中等温度(18.6℃)、中等DMC(35.6)、中等DC(220.7)
   • 物理解释：典型的季节性干燥条件
   • 实际意义：常规防火措施适用

3. 低风险集群（集群2） ：

   • 特征：低温(12.7℃)、低DMC(15.2)、低DC(110.5)
   • 物理解释：湿润凉爽条件
   • 实际意义：火灾风险最低，可减少防火资源投入

集群特征雷达图对比 ：

from math import pi

# 准备数据
categories = features
N = len(categories)
angles = [n / float(N) * 2 * pi for n in range(N)]
angles += angles[:1]

fig = plt.figure(figsize=(8, 8))
ax = fig.add_subplot(111, polar=True)
ax.set_theta_offset(pi / 2)
ax.set_theta_direction(-1)

plt.xticks(angles[:-1], categories)
ax.set_rlabel_position(0)
plt.yticks([0.5, 1, 1.5], ["0.5", "1.0", "1.5"], color="grey", size=7)
plt.ylim(0, 1.5)

for i in range(optimal_k):
    values = centers[i]
    values = (values - values.min()) / (values.max() - values.min())  # 归一化
    values = np.append(values, values[0])
    ax.plot(angles, values, linewidth=2, linestyle='solid', label=f'集群 {i}')
    ax.fill(angles, values, alpha=0.1)

plt.title('聚类中心特征对比雷达图', y=1.1)
plt.legend(loc='upper right', bbox_to_anchor=(1.3, 1.1))
plt.show()

6. 项目扩展与实用建议

完成基础分析后，可以考虑以下扩展方向：

1. 动态聚类过程可视化 ：

   from sklearn.datasets import make_blobs
   from IPython.display import clear_output
   import time
   
   # 生成模拟数据
   X, _ = make_blobs(n_samples=300, centers=3, random_state=42)
   
   # 初始化
   kmeans = KMeans(n_clusters=3, init='random', max_iter=1, random_state=42, n_init=1)
   
   plt.figure(figsize=(10, 8))
   for i in range(1, 10):
       kmeans.fit(X)
       clear_output(wait=True)
       
       plt.scatter(X[:, 0], X[:, 1], c=kmeans.labels_, cmap='viridis', s=50)
       plt.scatter(kmeans.cluster_centers_[:, 0], kmeans.cluster_centers_[:, 1], 
                   c='red', marker='X', s=200)
       plt.title(f'迭代次数: {i}')
       plt.show()
       time.sleep(1)
   

2. 结合时间序列分析 ：

   • 按月/季分析集群分布变化
   • 预测高风险时段

3. 部署为Web应用 ：

   from flask import Flask, render_template
   import plotly.express as px
   
   app = Flask(__name__)
   
   @app.route('/')
   def visualize():
       fig = px.scatter_3d(data, x='FFMC', y='DMC', z='Temp',
                          color='Cluster', opacity=0.7)
       return render_template('viz.html', plot=fig.to_html())
   
   if __name__ == '__main__':
       app.run(debug=True)
   

注意：在实际项目中，建议使用更高效的聚类算法如MiniBatchKMeans处理大规模数据，并考虑使用Plotly或Bokeh等交互性更强的可视化库。