# 02.2 空间推理

## 核心问题

> 机器如何理解和处理空间关系？
> 图算法在空间分析中有哪些应用？
> 如何进行连通性分析和路径优化？

---

## 概念讲解

### 空间关系类型

```
空间关系分类

┌─────────────────────────────────────────────────────────────┐
│                                                             │
│  1. 拓扑关系       │
│     - 相邻 (Adjacent): A与B共享边界                          │
│     - 包含 (Contains): A完全包含B                           │
│     - 重叠 (Overlaps): A与B部分重叠                         │
│     - 相离 (Disjoint): A与B不相交                           │
│                                                             │
│  2. 距离关系         │
│     - 欧氏距离: 直线距离                                     │
│     - 曼哈顿距离: 城市街区距离                               │
│     - 阻力距离: 穿越不同地形的代价                           │
│     - 时间距离: 行驶时间成本                                 │
│                                                             │
│  3. 方向关系       │
│     - 绝对方向: 北、南、东、西                               │
│     - 相对方向: 前、后、左、右                               │
│     - 方位角: 0-360度的精确方向                             │
│                                                             │
│  4. 模式关系         │
│     - 聚集: 要素密集分布                                     │
│     - 离散: 要素分散分布                                     │
│     - 随机: 要素随机分布                                     │
│     - 规则: 要素有规律分布                                  │
│                                                             │
└─────────────────────────────────────────────────────────────┘
```

### 图表示与空间推理

空间问题常可转换为图问题：

```
空间 → 图的转换

空间场景          图表示
───────────      ───────
源地A ──廊道──→ 源地B      节点A ──边──→ 节点B
   │                     │
   └──廊道──→ 源地C       └──边──→ 节点C
```

**空间问题的图抽象**：

| 空间问题 | 图表示 | 算法 |
|---------|--------|------|
| 最短路径 | 节点=位置，边=路径 | Dijkstra, A* |
| 连通性分析 | 节点=斑块，边=廊道 | BFS, DFS, 并查集 |
| 设施选址 | 节点=候选点，边=需求 | p-median, p-center |
| 覆盖问题 | 节点=服务点，边=覆盖范围 | 最大覆盖 |
| 网络流 | 节点=源/汇，边=管道 | 最大流最小割 |

---

## 设计原理

### 连通性分析

连通性是生态网络分析的核心：

```python
class ConnectivityAnalyzer:
    """
    连通性分析器

    核心：使用图算法分析空间连通性
    """

    def __init__(self, resistance_surface):
        """
        Args:
            resistance_surface: 阻力面栅格
        """
        self.resistance = resistance_surface
        self.graph = None

    def build_graph(self):
        """将阻力面转换为图"""
        import networkx as nx

        # 创建图
        self.graph = nx.Graph()

        rows, cols = self.resistance.shape

        # 添加节点和边
        for i in range(rows):
            for j in range(cols):
                node_id = i * cols + j

                # 添加节点
                self.graph.add_node(node_id, pos=(i, j))

                # 添加边（8邻域）
                for di in [-1, 0, 1]:
                    for dj in [-1, 0, 1]:
                        if di == 0 and dj == 0:
                            continue

                        ni, nj = i + di, j + dj
                        if 0 <= ni < rows and 0 <= nj < cols:
                            neighbor_id = ni * cols + nj

                            # 边权重 = 平均阻力
                            weight = (
                                self.resistance[i, j] +
                                self.resistance[ni, nj]
                            ) / 2

                            self.graph.add_edge(
                                node_id, neighbor_id,
                                weight=weight
                            )

        return self.graph

    def least_cost_path(self, source, target):
        """计算最小阻力路径"""
        if self.graph is None:
            self.build_graph()

        # Dijkstra算法
        path = nx.shortest_path(
            self.graph,
            source=source,
            target=target,
            weight='weight'
        )

        return path

    def connectivity_metrics(self, sources):
        """
        计算连通性指标

        Args:
            sources: 源地节点列表

        Returns:
            连通性指标字典
        """
        if self.graph is None:
            self.build_graph()

        metrics = {}

        # 1. 整体连通性 (图的连通分量数)
        components = list(nx.connected_components(
            self.graph.subgraph(sources)
        ))
        metrics['n_components'] = len(components)

        # 2. 最大连通分量大小
        if components:
            metrics['largest_component'] = max(len(c) for c in components)
        else:
            metrics['largest_component'] = 0

        # 3. 平均最短路径长度
        if len(sources) > 1:
            subgraph = self.graph.subgraph(sources)
            if nx.is_connected(subgraph):
                metrics['avg_path_length'] = nx.average_shortest_path_length(
                    subgraph, weight='weight'
                )
            else:
                metrics['avg_path_length'] = float('inf')

        # 4. 网络密度
        n = len(sources)
        if n > 1:
            max_edges = n * (n - 1) / 2
            actual_edges = self.graph.subgraph(sources).number_of_edges()
            metrics['density'] = actual_edges / max_edges
        else:
            metrics['density'] = 0

        return metrics
```

### 最短路径算法

空间分析中最常用的图算法：

```python
"""
最短路径算法比较
"""
import heapq
from typing import Dict, List, Tuple, Set

class ShortestPathAlgorithms:
    """最短路径算法集合"""

    def __init__(self, graph: Dict):
        """
        Args:
            graph: {node: {neighbor: weight, ...}, ...}
        """
        self.graph = graph

    def dijkstra(self, start: str, goal: str = None) -> Tuple[Dict, Dict]:
        """
        Dijkstra算法：经典最短路径

        适合：非负权重图
        复杂度：O((V+E)logV)
        """
        # 优先队列：(距离, 节点)
        pq = [(0, start)]
        visited = set()
        distances = {start: 0}
        parents = {start: None}

        while pq:
            current_dist, current = heapq.heappop(pq)

            if current in visited:
                continue
            visited.add(current)

            if current == goal:
                break

            for neighbor, weight in self.graph.get(current, {}).items():
                if neighbor in visited:
                    continue

                new_dist = current_dist + weight

                if new_dist < distances.get(neighbor, float('inf')):
                    distances[neighbor] = new_dist
                    parents[neighbor] = current
                    heapq.heappush(pq, (new_dist, neighbor))

        return distances, parents

    def reconstruct_path(self, parents: Dict, start: str, goal: str) -> List:
        """从parents字典重建路径"""
        path = []
        current = goal

        while current is not None:
            path.append(current)
            current = parents.get(current)

        path.reverse()

        if path[0] == start:
            return path
        return []

    def a_star(self, start: str, goal: str,
               heuristic: callable) -> Tuple[Dict, Dict]:
        """
        A*算法：带启发式的最短路径

        适合：有目标节点的图，有可用启发式
        复杂度：O(b^d) 实际通常比Dijkstra快
        """
        def h(node):
            return heuristic(node, goal)

        # f(n) = g(n) + h(n)
        pq = [(h(start), 0, start)]
        visited = set()
        g_score = {start: 0}  # 实际距离
        parents = {start: None}

        while pq:
            f, g, current = heapq.heappop(pq)

            if current in visited:
                continue
            visited.add(current)

            if current == goal:
                break

            for neighbor, weight in self.graph.get(current, {}).items():
                if neighbor in visited:
                    continue

                tentative_g = g + weight

                if tentative_g < g_score.get(neighbor, float('inf')):
                    g_score[neighbor] = tentative_g
                    f_score = tentative_g + h(neighbor)
                    parents[neighbor] = current
                    heapq.heappush(pq, (f_score, tentative_g, neighbor))

        return g_score, parents

# 空间启发式函数
def euclidean_heuristic(node_pos: Tuple, goal_pos: Tuple) -> float:
    """欧氏距离启发式"""
    import math
    return math.sqrt(
        (node_pos[0] - goal_pos[0])**2 +
        (node_pos[1] - goal_pos[1])**2
    )

def manhattan_heuristic(node_pos: Tuple, goal_pos: Tuple) -> float:
    """曼哈顿距离启发式（适合网格）"""
    return abs(node_pos[0] - goal_pos[0]) + abs(node_pos[1] - goal_pos[1])
```

---

## 代码示例

### 生态廊道识别

```python
"""
基于空间推理的生态廊道识别
"""
import numpy as np
from typing import List, Tuple
import heapq

def extract_corridors_mcr(resistance_surface: np.ndarray,
                         sources: List[Tuple[int, int]]) -> List[dict]:
    """
    使用最小累积阻力(MCR)方法提取生态廊道

    Args:
        resistance_surface: 阻力面栅格
        sources: 源地坐标列表 [(row, col), ...]

    Returns:
        廊道列表
    """
    rows, cols = resistance_surface.shape

    # 计算成本距离
    cost_distance = compute_cost_distance(resistance_surface, sources)

    # 提取廊道（低阻力通道）
    corridors = []

    for i, source1 in enumerate(sources):
        for source2 in sources[i+1:]:
            # 找到两源之间的最低阻力路径
            path = extract_lowest_resistance_path(
                cost_distance, resistance_surface, source1, source2
            )

            if path:
                corridors.append({
                    'source_a': source1,
                    'source_b': source2,
                    'path': path,
                    'cost': sum(resistance_surface[p] for p in path)
                })

    return corridors

def compute_cost_distance(resistance: np.ndarray,
                         sources: List[Tuple[int, int]]) -> np.ndarray:
    """
    计算成本距离（到最近源地的累积阻力）

    使用Dijkstra算法的变种
    """
    rows, cols = resistance.shape
    cost = np.full((rows, cols), np.inf)

    # 优先队列：(累积成本, row, col)
    pq = []

    # 初始化源地
    for source_row, source_col in sources:
        cost[source_row, source_col] = 0
        heapq.heappush(pq, (0, source_row, source_col))

    # 8方向
    directions = [(-1, 0), (1, 0), (0, -1), (0, 1),
                  (-1, -1), (-1, 1), (1, -1), (1, 1)]

    visited = np.zeros((rows, cols), dtype=bool)

    while pq:
        current_cost, row, col = heapq.heappop(pq)

        if visited[row, col]:
            continue
        visited[row, col] = True

        for dr, dc in directions:
            nr, nc = row + dr, col + dc

            if 0 <= nr < rows and 0 <= nc < cols:
                # 计算移动成本
                if dr != 0 and dc != 0:  # 对角移动
                    move_cost = resistance[nr, nc] * 1.414
                else:
                    move_cost = resistance[nr, nc]

                new_cost = current_cost + move_cost

                if new_cost < cost[nr, nc]:
                    cost[nr, nc] = new_cost
                    heapq.heappush(pq, (new_cost, nr, nc))

    return cost

def extract_lowest_resistance_path(cost_distance: np.ndarray,
                                   resistance: np.ndarray,
                                   start: Tuple[int, int],
                                   end: Tuple[int, int]) -> List[Tuple[int, int]]:
    """
    从成本距离表面提取最低阻力路径
    """
    path = [end]
    current = end

    while current != start:
        row, col = current
        best_neighbor = None
        best_cost = cost_distance[current]

        # 检查邻域
        for dr in [-1, 0, 1]:
            for dc in [-1, 0, 1]:
                if dr == 0 and dc == 0:
                    continue

                nr, nc = row + dr, col + dc
                if (0 <= nr < cost_distance.shape[0] and
                    0 <= nc < cost_distance.shape[1]):
                    if cost_distance[nr, nc] < best_cost:
                        best_cost = cost_distance[nr, nc]
                        best_neighbor = (nr, nc)

        if best_neighbor is None:
            break

        path.append(best_neighbor)
        current = best_neighbor

    path.reverse()
    return path if path[0] == start else []
```

---

## 案例分析

### ENAgent中的廊道识别

ENAgent使用空间推理提取生态廊道：

```python
class ENAgentCorridorExtractor:
    """ENAgent的廊道提取模块"""

    def extract_corridors(self, mcr_surface, sources, width_threshold=500):
        """
        基于MCR表面提取廊道

        Args:
            mcr_surface: 最小累积阻力表面
            sources: 源地列表
            width_threshold: 廊道最小宽度

        Returns:
            廊道字典
        """
        corridors = []

        # 对每对源地提取路径
        for i in range(len(sources)):
            for j in range(i + 1, len(sources)):
                path = self._extract_path_between_sources(
                    mcr_surface, sources[i], sources[j]
                )

                if path:
                    # 分析廊道宽度
                    width = self._calculate_corridor_width(
                        mcr_surface, path
                    )

                    if width >= width_threshold:
                        corridors.append({
                            'from': sources[i]['id'],
                            'to': sources[j]['id'],
                            'path': path,
                            'width': width,
                            'quality': self._assess_quality(
                                mcr_surface, path
                            )
                        })

        return corridors
```

---

## 反思与延伸

### 思考问题

1. **算法选择**：什么时候用Dijkstra，什么时候用A*？

2. **空间尺度**：空间推理如何处理多尺度问题？

3. **计算效率**：大规模空间数据的图算法如何优化？

4. **动态变化**：空间环境变化时，如何高效更新推理结果？

### 延伸阅读

- **"Network Flows"** (Ahuja, Magnanti, Orlin) - 网络流理论
- **"Geometric Algorithms"** - 几何算法
- NetworkX文档 - Python图算法库

---

## 关键要点

1. **空间关系有四类**：拓扑、距离、方向、模式
2. **图算法是空间推理的核心工具**
3. **连通性分析**使用图的结构特性
4. **最短路径**有多个算法变种，各有适用场景
5. **MCR分析**本质是图上的最短路径问题