### Python 实现 ISODATA 聚类算法 以下是基于 Python 的 ISODATA 聚类算法实现，适用于遥感图像数据处理。此代码依赖于 NumPy 和 SciPy 库来完成矩阵运算和距离计算。 #### 1. 导入必要的库 首先导入所需的第三方库： ```python import numpy as np from scipy.spatial.distance import cdist ``` #### 2. 定义 ISODATA 类 定义一个名为 `ISODATA` 的类，封装算法的主要逻辑： ```python class ISODATA: def __init__(self, data, initial_clusters=8, max_iterations=100, min_samples_per_cluster=5, split_threshold=10, merge_threshold=5): self.data = data self.initial_clusters = initial_clusters self.max_iterations = max_iterations self.min_samples_per_cluster = min_samples_per_cluster self.split_threshold = split_threshold self.merge_threshold = merge_threshold # 初始化簇中心 indices = np.random.choice(len(data), size=self.initial_clusters, replace=False) self.centers = data[indices] def assign_points(self): distances = cdist(self.data, self.centers, 'euclidean') labels = np.argmin(distances, axis=1) return labels def update_centers(self, labels): new_centers = [] counts = [] for i in range(len(self.centers)): cluster_data = self.data[labels == i] if len(cluster_data) > 0: center = np.mean(cluster_data, axis=0) new_centers.append(center) counts.append(len(cluster_data)) else: counts.append(0) self.centers = np.array(new_centers) return counts def split_clusters(self, counts): to_split = [ idx for idx, count in enumerate(counts) if count >= self.split_threshold and len(self.centers[idx]) > 1 ] new_centers = list(self.centors) for idx in to_split: original_center = self.centers[idx] offset = np.std(original_center) * np.random.randn(len(original_center)) / 2 new_centers.remove(original_center) new_centers.extend([original_center + offset, original_center - offset]) self.centers = np.array(new_centers) def merge_clusters(self): pairwise_distances = cdist(self.centers, self.centers, 'euclidean') np.fill_diagonal(pairwise_distances, np.inf) while np.any(pairwise_distances < self.merge_threshold): closest_pair = np.unravel_index(np.argmin(pairwise_distances), pairwise_distances.shape) merged_center = (self.centers[closest_pair[0]] + self.centers[closest_pair[1]]) / 2 mask = np.ones(len(self.centers), dtype=bool) mask[[closest_pair[0], closest_pair[1]]] = False self.centers = np.vstack((self.centers[mask], merged_center)) pairwise_distances = cdist(self.centers, self.centers, 'euclidean') np.fill_diagonal(pairwise_distances, np.inf) def fit(self): iteration = 0 while iteration < self.max_iterations: labels = self.assign_points() counts = self.update_centers(labels) self.split_clusters(counts) self.merge_clusters() iteration += 1 return labels ``` #### 3. 使用示例 加载遥感数据并调用 ISODATA 算法进行聚类分析： ```python # 假设我们有一个形状为 (n_samples, n_features) 的遥感数据数组 data = np.random.rand(1000, 3) # 示例数据 model = ISODATA(data=data, initial_clusters=5, max_iterations=50, min_samples_per_cluster=10, split_threshold=20, merge_threshold=5) cluster_labels = model.fit() # 运行模型得到标签 print("Cluster Labels:", cluster_labels) ``` --- ### 关键点说明 - **初始化阶段**：通过随机采样初始样本作为簇中心[^1]。 - **分配阶段**：利用欧几里得距离将每个样本指派到最近的簇中心[^2]。 - **更新阶段**：重新计算各簇的新中心位置，并统计每簇中的样本数量[^3]。 - **分裂与合并机制**：当某个簇满足特定条件时会触发分裂；反之，如果两个簇之间的距离小于阈值，则执行合并操作[^2]。 --- ###