TowardsDataScience 2023 博客中文翻译（七十八）

原文：TowardsDataScience

协议：CC BY-NC-SA 4.0

完整的 Python 缓存指南

原文：towardsdatascience.com/complete-guide-to-caching-in-python-b4e37a4bcebf

缓存如何工作及缓存函数的方法

Kay Jan Wong

·发布于 Towards Data Science ·阅读时间 7 分钟·2023 年 12 月 1 日

–

照片由 Nana Smirnova 提供，来源于 Unsplash

当重复调用函数并使用相同参数时，会导致计算被重复执行。记忆化 在这种情况下非常有用，因为函数调用的结果可以被‘保存’以供将来使用。这将节省时间并优化代码，因为代码变得更少计算密集。缓存 是一个更通用的术语，用于指代存储任何数据。

本文将探讨不同的缓存策略、缓存考虑因素以及如何启用和实现不同类型的缓存（使用 Python 包和你的实现）！

目录

• 缓存类型

• 缓存考虑因素

• LRU 缓存

• LFU 缓存

• FIFO / LIFO 缓存

缓存类型

根据你的需求，有几种缓存策略，例如：

• 最少最近使用（LRU）：移除最少最近使用的数据，是最常见的缓存类型

• 最少频繁使用（LFU）：移除最少频繁使用的数据

• 先进先出（FIFO）：移除最旧的数据

• 后进先出（LIFO）：移除最新的数据

• 最近使用（MRU）：移除最近使用的数据

• 随机替换（RR）：移除随机选择的数据

缓存考虑因素

在应用程序中使用缓存时，您应该考虑缓存的内存占用，因为它存储了额外的信息。如果在不同的实现之间进行选择，就架构和数据结构而言，有几个时间考虑因素，例如：

• 访问时间：对于已经计算过的参数，结果应在O(1)时间内快速访问

• 插入时间：对于新参数，数据应尽可能在O(1)时间内插入缓存（根据实现，有些可能需要O(n)时间，选择时请谨慎！）

• 删除时间：当缓存满时，需要根据缓存策略删除数据。删除涉及识别要删除的数据，然后从内存中移除它们。

1 — 最少使用（LRU）缓存

通过为 Python 函数添加装饰器来实现

图 1：LRU 实现（图片由作者提供）

工作原理：它使用字典和双向链表。在字典中，键值对是提供的参数和双向链表中的条目。这使得根据提供的参数快速引用结果（O(1)访问时间）。由于参数可以作为字典键存储，因此它们必须是可哈希的。

在双向链表中，函数结果存储在条目中。条目按其最近性排序，并可以引用其直接前一个和下一个条目。这允许条目轻松插入或重新排序，并在缓存满时快速识别和删除最久未用的条目。

from functools import lru_cache

@lru_cache
def fibonacci(n: int) -> int:
    return n if n < 2 else fibonacci(n - 1) + fibonacci(n - 2)

高级用法

通过实现以下功能来增强 LRU 缓存：

1. 使用maxsize参数设置最大缓存大小

2. 使用maxsize=None参数或functools.cache设置无限缓存大小。这意味着数据将永远不会被删除，这可能在内存占用的代价下带来时间改进。

3. 使用.cache_info()检索缓存命中和未命中的信息。缓存统计信息使我们能够评估缓存的访问效率。

4. 使用cachetools.TTLCache添加过期时间

5. 将缓存与 CPU 内存使用挂钩而不是缓存大小

# 1\. Maximum cache size
@lru_cache(maxsize=1024)
def fibonacci(n: int) -> int:
    return n if n < 2 else fibonacci(n - 1) + fibonacci(n - 2)

# 2\. Unlimited cache size
@lru_cache(maxsize=None)
def fibonacci(n: int) -> int:
    return n if n < 2 else fibonacci(n - 1) + fibonacci(n - 2)

from functools import cache

@cache
def fibonacci(n: int) -> int:
    return n if n < 2 else fibonacci(n - 1) + fibonacci(n - 2)

# 3\. Retrieve caching information
fibonacci.cache_info()

# 4\. Add expiration time (requires pip install cachetools)
from cachetools import cached, TTLCache

@cached(cache=TTLCache(maxsize=100, ttl=60*60*5), info=True)
def fibonacci(n: int) -> int:
    return n if n < 2 else fibonacci(n - 1) + fibonacci(n - 2)

注意：其他 Python 包实现了缓存，例如fastcache，但它们没有functools.lru_cache那么受欢迎或常用。

2 — 最少频繁使用（LFU）缓存

通过维护一个缓存字典和一个频率字典来实现

在线可以找到许多 LFU 缓存的实现，包括cachedtools，其使用方法类似于上一节中的 LRU 缓存。LFU 缓存可以使用哈希映射、单向链表或双向链表实现。

我发现维护两个字典是最优的方式，考虑到访问、插入和删除时间。通过使用哈希可以进一步提高内存使用效率。

图 2: LFU 实现（作者图片）

工作原理：它使用 2 个字典。在缓存字典中，键是提供的参数，值是一个包含函数结果和频率的元组。这允许快速检索函数结果（O(1)访问时间！）和频率（用于访问频率字典）。

在频率字典中，键是频率，值是提供的参数列表。存储提供的参数列表允许快速识别最不常用的参数，并从两个字典中逐出该参数及其结果。可以使用 Deque 代替列表，以便更快的访问、插入和删除时间。

警告：在 LFU 缓存中，它对最近的条目有偏见，因为较新的条目可能没有现有条目的频率高——使得即使访问频率较低，也很难逐出较旧的条目。

from collections import defaultdict, deque
from typing import Any, Deque, Dict, Tuple

class LFUCache:
    def __init__(self, maxsize: int | None = None):
        """
        Args:
            maxsize (int | None): Capacity of cache size, defaults to None
        """
        if maxsize and maxsize < 0:
            maxsize = 0
        self.maxsize = maxsize
        self.cache_dict: Dict[Any, Tuple[Any, int]] = {}
        self.freq_dict: Dict[int, Deque[Any]] = defaultdict(lambda: deque([]))
        self.hits, self.misses, self.curr_size = 0, 0, 0

    def cache_info(self) -> Dict[str, int | None]:
        """Report cache statistics"""
        return dict(
            hits=self.hits,
            misses=self.misses,
            maxsize=self.maxsize,
            currsize=self.curr_size,
        )

    def cache_clear(self) -> None:
        """Clear the cache and cache statistics"""
        self.cache_dict = {}
        self.freq_dict = defaultdict(lambda: deque([]))
        self.hits, self.misses, self.curr_size = 0, 0, 0

    def update(self, key: Any, value: Any) -> None:
        """Update frequency of key in cache and frequency dictionary.
        Removes key in frequency dictionary if necessary.

        Args:
            key (Any): Argument to function
            value (Any): Result of function
        """
        _, freq = self.cache_dict[key]
        self.cache_dict[key] = (value, freq + 1)
        self.freq_dict[freq].remove(key)
        self.freq_dict[freq + 1].append(key)
        if not len(self.freq_dict[freq]):
            del self.freq_dict[freq]

    def get(self, key: Any) -> Any:
        """Get value by key. Updates the hits and misses statistics.

        Args:
            key (Any): Argument to function

        Returns:
            (Any)
        """
        if key in self.cache_dict:
            self.hits += 1
            value, _ = self.cache_dict[key]
            self.update(key, value)
            return value
        self.misses += 1
        raise KeyError(f"{key} does not exist in cache.")

    def put(self, key: Any, value: Any) -> None:
        """Put value by key into cache and frequency dictionary.
        Check the capacity of the cache and delete the key-value if necessary.

        Args:
            key (Any): Argument to function
            value (Any): Result of function
        """
        if key in self.cache_dict:
            self.update(key, value)
        else:
            self.cache_dict[key] = (value, 1)
            self.freq_dict[1].append(key)
            self.curr_size += 1
            if self.maxsize is not None and self.curr_size > self.maxsize:
                remove_key = self.freq_dict[min(self.freq_dict)].popleft()
                del self.cache_dict[remove_key]
                self.curr_size -= 1

上面的代码片段改编自这里，并做了一些调整。

要将 LFU 缓存作为装饰器使用，我们可以封装LFUCache类，并类似于functools.lru_cache来使用。可以按如下方式进行：

from functools import wraps
from typing import Callable

def lfu_cache(maxsize: int | None = None) -> Any:
    cache = LFUCache(maxsize)

    def decorator(func: Callable[..., Any]) -> Any:

        @wraps(func)
        def wrapper(*args, **kwargs) -> Callable[..., Any]:
            key = hash(*args, **kwargs)
            try:
                result = cache.get(key)
                return result
            except KeyError:
                result = func(*args, **kwargs)
                cache.put(key, result)
            return result

        wrapper.cache = cache
        wrapper.cache_info = cache.cache_info
        wrapper.cache_clear = cache.cache_clear

        return wrapper
    return decorator

# Usage
@lfu_cache(maxsize=1024)
def fibonacci(n: int) -> int:
    return n if n < 2 else fibonacci(n - 1) + fibonacci(n - 2)

fibonacci.cache_info()

3 — 先进先出 (FIFO) / 后进先出 (LIFO) 缓存

通过维护一个有序字典实现

图 3: FIFO / LIFO 实现（作者图片）

工作原理：它使用一个有序字典，其中键值对是提供的参数和函数结果。字典按照提供的参数首次调用的时间进行排序。

这允许根据提供的参数快速引用结果（O(1)访问时间），条目可以根据缓存策略从前面或后面移除（O(1)插入和删除时间）。实现可以使用cachedtools，使用方式类似于上一节中的 LRU 缓存。为了说明这一点：

# Requires pip install cachetools
from cachetools import cached, FIFOCache

@cached(cache=FIFOCache(maxsize=100), info=True)
def fibonacci(n: int) -> int:
    return n if n < 2 else fibonacci(n - 1) + fibonacci(n - 2)

# Retrieve caching information
fibonacci.cache_info()

希望你对缓存、缓存策略的类型及其考虑因素，以及使用 Python 库或你自己的实现有了更多了解！

相关链接

• LRU 文档: docs.python.org/3/library/functools.html#functools.lru_cache

• LFU 文章: www.tutorialspoint.com/lfu-cache-in-python

• cachetools 官方 GitHub: github.com/tkem/cachetools

完整实现一个用于图像识别的迷你 VGG 网络

原文：towardsdatascience.com/complete-implementation-of-a-mini-vgg-network-for-image-recognition-849299480356

照片由 Guillaume de Germain 提供，来源于 Unsplash

用于更高效图像识别的深度卷积神经网络

Rashida Nasrin Sucky

·发表于 Towards Data Science ·7 分钟阅读·2023 年 2 月 27 日

–

VGG 网络是最流行的图像识别技术之一的基础。它值得学习，因为它打开了许多新的领域。要理解 VGGNet，你需要了解卷积神经网络（CNN）。如果你不熟悉 CNN 架构，请先阅读本教程：

## 卷积神经网络：对层的良好理解及图像分类示例

内容丰富

towardsdatascience.com

在本文中，我们将只关注 VGGNet 的实现部分。因此，我们将在这里快速进行。

关于 VGG 网络

VGGNet 是一种可以更成功地提取特征的卷积神经网络（CNN）。在 VGGNet 中，我们堆叠了多个卷积层。VGGNet 可以是浅层的，也可以是深层的。在浅层 VGGNet 中，通常仅添加两组四个卷积层，如我们将很快看到的那样。而在深层 VGGNet 中，可以添加超过四个卷积层。两个常用的深层 VGGNet 是 VGG16，它总共使用 16 层，以及 VGG19，它总共使用 19 层。我们可以添加一个批量归一化层，也可以不添加。但在本教程中我将使用它。

你可以通过这个链接阅读更多关于架构的信息：

[## VGG 非常深的卷积网络 (VGGNet) - 你需要知道的 - viso.ai

我们使用 cookie 来提升您的浏览体验、提供个性化广告或内容，并分析我们的流量。通过…

viso.ai](https://viso.ai/deep-learning/vgg-very-deep-convolutional-networks/?source=post_page-----849299480356--------------------------------)

今天我们将研究 mini VGGNet。因此，它会更简单、更容易运行，但对于许多用例仍然很强大。

miniVGGNet 的一个重要特点是它使用所有的 3x3 滤波器。这就是它能够如此出色地泛化的原因。让我们开始在 Keras 和 TensorFlow 中构建一个 mini VGGNet。

我使用了 Google Colaboratory 笔记本，并启用了 GPU。否则，训练速度非常慢。

Mini VGG 网络开发、训练和评估

现在开始工作吧。我们将对其进行一些实验，以展示如何使用它。

这些是必要的导入：

import tensorflow as tf
from keras.models import Sequential
from keras.layers.normalization import BatchNormalization
from keras.layers.convolutional import Conv2D
from keras.layers.convolutional import MaxPooling2D
from keras.layers.core import Activation
from keras.layers.core import Flatten
from keras.layers.core import Dropout
from keras.layers.core import Dense
from keras import backend as K
from sklearn.preprocessing import LabelBinarizer
from sklearn.metrics import classification_report
from keras.optimizers import SGD
from keras.datasets import cifar10
import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline

这些导入文件真不少！

我们将使用cifar-10 数据集，这是一个在 TensorFlow 库中提供的公共数据集。

我使用了两个不同的网络作为实验。第一个是流行的。我之所以说流行，是因为我在Kaggle和其他一些教程中找到了这个[架构]。

class MiniVGGNet:
 @staticmethod
 def build(width, height, depth, classes):
  # initialize the model along with the input shape to be
  # "channels last" and the channels dimension itself
  model = Sequential()
  inputShape = (height, width, depth)
  chanDim = -1

  if K.image_data_format() == "channels_first":
   inputShape = (depth, height, width)
   chanDim = 1

  # first CONV => Activation => CONV => Activation => POOL layer set
  model.add(Conv2D(32, (3, 3), padding="same",
   input_shape=inputShape))
  model.add(Activation("relu"))
  model.add(BatchNormalization(axis=chanDim))
  model.add(Conv2D(32, (3, 3), padding="same"))
  model.add(Activation("relu"))
  model.add(BatchNormalization(axis=chanDim))
  model.add(MaxPooling2D(pool_size=(2, 2)))
  model.add(Dropout(0.25))

  # second CONV => Activation => CONV => Activation => POOL layer set
  model.add(Conv2D(64, (3, 3), padding="same"))
  model.add(Activation("relu"))
  model.add(BatchNormalization(axis=chanDim))
  model.add(Conv2D(64, (3, 3), padding="same"))
  model.add(Activation("relu"))
  model.add(BatchNormalization(axis=chanDim))
  model.add(MaxPooling2D(pool_size=(2, 2)))
  model.add(Dropout(0.25))

  # Dense Layer
  model.add(Flatten())
  model.add(Dense(512))
  model.add(Activation("relu"))
  model.add(BatchNormalization())
  model.add(Dropout(0.5))

  # softmax classifier
  model.add(Dense(classes))
  model.add(Activation("softmax"))

  # return the constructed network architecture
  return model

让我们加载并准备我们的 cifar-10 数据集。

(x_train, y_train), (x_test, y_test) = cifar10.load_data()
x_train = x_train.astype("float") / 255.0 
x_test = x_test.astype("float") / 255.0 

cifar-10 数据集有 10 个标签。这些是 cifar-10 数据集中的标签：

labelNames = ["airplane", "automobile", "bird", "cat", "deer",
 "dog", "frog", "horse", "ship", "truck"]

使用 LabelBinarizer 将标签二值化：

lb = LabelBinarizer()
y_train = lb.fit_transform(y_train)
y_test = lb.transform(y_test)

在这里编译模型。评估指标是“准确率”，我们将运行 10 个周期。

optimizer = tf.keras.optimizers.legacy.SGD(learning_rate=0.01, decay=0.01/40, momentum=0.9,
                                           nesterov=True)
model = miniVGGNet.build(width = 32, height = 32, depth = 3, classes=10)
model.compile(loss='categorical_crossentropy', optimizer = optimizer,
              metrics=['accuracy'])
h = model.fit(x_train, y_train, validation_data=(x_test, y_test),
              batch_size = 64, epochs=10, verbose=1)

结果如下：

Epoch 1/10
782/782 [==============================] - 424s 539ms/step - loss: 1.6196 - accuracy: 0.4592 - val_loss: 1.4083 - val_accuracy: 0.5159
Epoch 2/10
782/782 [==============================] - 430s 550ms/step - loss: 1.1437 - accuracy: 0.6039 - val_loss: 1.0213 - val_accuracy: 0.6505
Epoch 3/10
782/782 [==============================] - 430s 550ms/step - loss: 0.9634 - accuracy: 0.6618 - val_loss: 0.8495 - val_accuracy: 0.7013
Epoch 4/10
782/782 [==============================] - 427s 546ms/step - loss: 0.8532 - accuracy: 0.6998 - val_loss: 0.7881 - val_accuracy: 0.7215
Epoch 5/10
782/782 [==============================] - 425s 543ms/step - loss: 0.7773 - accuracy: 0.7280 - val_loss: 0.8064 - val_accuracy: 0.7228
Epoch 6/10
782/782 [==============================] - 421s 538ms/step - loss: 0.7240 - accuracy: 0.7451 - val_loss: 0.6757 - val_accuracy: 0.7619
Epoch 7/10
782/782 [==============================] - 420s 537ms/step - loss: 0.6843 - accuracy: 0.7579 - val_loss: 0.6564 - val_accuracy: 0.7715
Epoch 8/10
782/782 [==============================] - 420s 537ms/step - loss: 0.6405 - accuracy: 0.7743 - val_loss: 0.6833 - val_accuracy: 0.7706
Epoch 9/10
782/782 [==============================] - 422s 540ms/step - loss: 0.6114 - accuracy: 0.7828 - val_loss: 0.6188 - val_accuracy: 0.7848
Epoch 10/10
782/782 [==============================] - 421s 538ms/step - loss: 0.5799 - accuracy: 0.7946 - val_loss: 0.6166 - val_accuracy: 0.7898

经过 10 个周期后，训练数据的准确率为 79.46%，验证数据的准确率为 78.98%。

考虑到这一点，我想对这个网络做一些更改，看看结果如何。让我们重新定义上面的网络。我在整个网络中使用了 64 个滤波器，密集层中有 256 个神经元，最后一个 dropout 层中有 40%的 dropout。

这是新的 mini VGG 网络：

class miniVGGNet:
  @staticmethod 

  def build(width, height, depth, classes):
    model = Sequential()
    inputShape = (height, width, depth)
    chanDim = -1 

    if K.image_data_format() == "channels_first":
      inputShape = (depth, height, width)
      chanDim = 1

    # first Conv => Activation => Conv => Activation => Pool layer set
    model.add(Conv2D(64, (3, 3), padding="same",
   input_shape=inputShape))
    model.add(Activation("relu"))
    model.add(BatchNormalization(axis=chanDim))
    model.add(Conv2D(64, (3, 3), padding="same"))
    model.add(Activation("relu"))
    model.add(BatchNormalization(axis=chanDim))
    model.add(MaxPooling2D(pool_size=(2, 2)))
    model.add(Dropout(0.25))

  # second Conv => Activation => Conv => Activation => Pool layer set
    model.add(Conv2D(64, (3, 3), padding="same"))
    model.add(Activation("relu"))
    model.add(BatchNormalization(axis=chanDim))
    model.add(Conv2D(64, (3, 3), padding="same"))
    model.add(Activation("relu"))
    model.add(BatchNormalization(axis=chanDim))
    model.add(MaxPooling2D(pool_size=(2, 2)))
    model.add(Dropout(0.25))

  # Dense Layer
    model.add(Flatten())
    model.add(Dense(300))
    model.add(Activation("relu"))
    model.add(BatchNormalization())
    model.add(Dropout(0.4))
    model.add(Dense(classes))
    model.add(Activation("softmax"))
    return model

我们将使用相同的优化参数和运行模型。但这里我使用了 20 个周期。

optimizer = tf.keras.optimizers.legacy.SGD(learning_rate=0.01, decay=0.01/40, momentum=0.9,
                                           nesterov=True)
model = miniVGGNet.build(width = 32, height = 32, depth = 3, classes=10)
model.compile(loss='categorical_crossentropy', optimizer = optimizer,
              metrics=['accuracy'])
h = model.fit(x_train, y_train, validation_data=(x_test, y_test),
              batch_size = 64, epochs=20, verbose=1)

结果如下：

Epoch 1/20
782/782 [==============================] - 22s 18ms/step - loss: 1.5210 - accuracy: 0.4697 - val_loss: 1.1626 - val_accuracy: 0.5854
Epoch 2/20
782/782 [==============================] - 14s 18ms/step - loss: 1.0706 - accuracy: 0.6219 - val_loss: 0.9913 - val_accuracy: 0.6586
Epoch 3/20
782/782 [==============================] - 14s 18ms/step - loss: 0.8947 - accuracy: 0.6826 - val_loss: 0.8697 - val_accuracy: 0.6941
Epoch 4/20
782/782 [==============================] - 14s 18ms/step - loss: 0.7926 - accuracy: 0.7208 - val_loss: 0.7649 - val_accuracy: 0.7294
Epoch 5/20
782/782 [==============================] - 14s 18ms/step - loss: 0.7192 - accuracy: 0.7470 - val_loss: 0.6937 - val_accuracy: 0.7593
Epoch 6/20
782/782 [==============================] - 13s 17ms/step - loss: 0.6641 - accuracy: 0.7640 - val_loss: 0.6899 - val_accuracy: 0.7639
Epoch 7/20
782/782 [==============================] - 13s 17ms/step - loss: 0.6141 - accuracy: 0.7805 - val_loss: 0.6589 - val_accuracy: 0.7742
Epoch 8/20
782/782 [==============================] - 13s 17ms/step - loss: 0.5774 - accuracy: 0.7960 - val_loss: 0.6565 - val_accuracy: 0.7734
Epoch 9/20
782/782 [==============================] - 14s 17ms/step - loss: 0.5430 - accuracy: 0.8077 - val_loss: 0.6092 - val_accuracy: 0.7921
Epoch 10/20
782/782 [==============================] - 14s 18ms/step - loss: 0.5145 - accuracy: 0.8177 - val_loss: 0.5904 - val_accuracy: 0.7944
Epoch 11/20
782/782 [==============================] - 13s 17ms/step - loss: 0.4922 - accuracy: 0.8256 - val_loss: 0.6041 - val_accuracy: 0.7975
Epoch 12/20
782/782 [==============================] - 14s 18ms/step - loss: 0.4614 - accuracy: 0.8381 - val_loss: 0.5889 - val_accuracy: 0.7981
Epoch 13/20
782/782 [==============================] - 14s 18ms/step - loss: 0.4358 - accuracy: 0.8457 - val_loss: 0.5590 - val_accuracy: 0.8120
Epoch 14/20
782/782 [==============================] - 13s 17ms/step - loss: 0.4186 - accuracy: 0.8508 - val_loss: 0.5555 - val_accuracy: 0.8092
Epoch 15/20
782/782 [==============================] - 13s 17ms/step - loss: 0.4019 - accuracy: 0.8582 - val_loss: 0.5739 - val_accuracy: 0.8108
Epoch 16/20
782/782 [==============================] - 14s 17ms/step - loss: 0.3804 - accuracy: 0.8658 - val_loss: 0.5577 - val_accuracy: 0.8136
Epoch 17/20
782/782 [==============================] - 13s 17ms/step - loss: 0.3687 - accuracy: 0.8672 - val_loss: 0.5544 - val_accuracy: 0.8170
Epoch 18/20
782/782 [==============================] - 13s 17ms/step - loss: 0.3541 - accuracy: 0.8744 - val_loss: 0.5435 - val_accuracy: 0.8199
Epoch 19/20
782/782 [==============================] - 13s 17ms/step - loss: 0.3438 - accuracy: 0.8758 - val_loss: 0.5533 - val_accuracy: 0.8167
Epoch 20/20
782/782 [==============================] - 13s 17ms/step - loss: 0.3292 - accuracy: 0.8845 - val_loss: 0.5491 - val_accuracy: 0.8199

如果你注意到，经过 10 个周期后，准确率略高于之前的网络，经过 20 个周期后准确率非常好。训练数据的准确率为 88.45%，验证数据的准确率为 81.99%。

在同一图表中展示训练和验证准确率以及训练和验证损失：

%matplotlib inline
plt.close('all')
plt.style.use("ggplot")
plt.figure(figsize=(8, 6))
plt.plot(np.arange(0, 20), h.history["loss"], label="train_loss")
plt.plot(np.arange(0, 20), h.history["val_loss"], label="val_loss")
plt.plot(np.arange(0, 20), h.history["accuracy"], label="train_acc")
plt.plot(np.arange(0, 20), h.history["val_accuracy"], label="val_acc")
plt.title("Training Loss and Accuracy")
plt.xlabel("No of Epochs")
plt.ylabel("Loss/Accuracy")
plt.legend()
plt.show()

作者提供的图片

训练损失平稳下降，验证损失也有一些波动地下降。

结论

请随意进行实验。根据项目尝试不同的参数，看看效果如何。我们稍后将研究一个深度网络。

随时在 Twitter 上关注我，并点赞我的 Facebook 页面。

更多阅读

OpenCV Python 中的简单边缘检测方法 [## OpenCV Python 中的简单边缘检测方法

高效使用 Canny 边缘检测

如何使用 OpenCV 进行阈值分割 [## 如何使用 OpenCV 进行阈值分割

简单的、Otsu 和自适应阈值实现示例

如何定义自定义层、激活函数和损失函数 [## 如何定义自定义层、激活函数和损失函数

逐步讲解和完整代码示例

如何定义自定义层、激活函数和损失函数 [## 逐步教程：在 TensorFlow 中开发多输出模型

包含完整代码

在 TensorFlow 中使用顺序 API 和函数 API 进行回归 [## TensorFlow 中使用顺序 API 和函数 API 进行回归

展示几种不同类型的模型结构

回归模型的 TensorFlow 实现

复利与指数分布

原文：towardsdatascience.com/compound-interest-and-the-exponential-distribution-402d5e7f446a

你的抵押贷款是无记忆的，指数分布也是如此

Rohit Pandey

·发表于 Towards Data Science ·阅读时长 8 分钟·2023 年 5 月 24 日

–

使用开源许可证通过 midjourney 生成的图像。

欧拉数、指数分布与复利

在某些时刻会发生许多有趣的事件。例如，公交车到达公交站、高速公路上的事故、足球比赛中的进球。模拟这些时刻事件的过程称为“点过程”。这些过程中的一个重要考虑因素是从一个事件到下一个事件所需的时间。例如，假如你刚刚错过了一辆公交车，你还要等多久才能等到下一辆公交车？这个时间是一个随机变量，随机变量的选择决定了点过程。这个随机变量的一个选择是一些不是很随机的东西（确定性数字只是随机数字的一个特殊情况）。例如，公交车按时刻表到达，每 10 分钟一次。这可能听起来像是最简单的点过程，但还有更简单的情况。当事件之间的时间遵循指数分布时（这个过程叫做泊松过程），它就是一个简单的过程。它被称为指数分布是有原因的。它与欧拉数e和复利息有关。在这篇文章中，我们将探讨它们之间的联系。

复利

假设你在银行存入 1 美元。年利率为x。在年末，你的余额将是*(1+x)。为了获得更多的钱，你要求银行每月支付利息而不是每年支付。由于利率为每年x*，因此你在一个月内赚取的利息将是x/12。你会立即将这笔利息再投资。因此，在第二个月，你的投资变为*(1+x/12)，这个数额再次按(1+x/12)增长，这意味着两个月后的金额是(1+x/12)²*。重复 12 个月，你年末的余额将是*(1+x/12)¹²*。利用二项式定理，年末的新余额是：

我们可以看到这比我们之前得到的*(1+x)更多。这是有道理的，因为我们在几个月内获得了利息，利息被再投资并进一步赚取了利息。但是为什么只停留在12* 个间隔呢？你希望尽可能频繁地进行复利。如果银行允许的话，每毫秒一次。我们将 12 个间隔推广到n个间隔，并使n非常大。在每个间隔后，我们的余额进一步增长*(1+x/n)^n*。在年末，我们将拥有的金额，

Eq (1)

使用二项式定理展开，

当n变得更大时，n-1、n-2 等几乎与n相同。因此，所有涉及n的项在分子和分母之间相互抵消（因为我们有n→∞），剩下的是：

如果我们对x求导数，则会得到B(x)本身。如果我们将x=1代入，会得到一个非常特殊的数字。你能猜到是什么吗？从前两项很明显，这个数字大于2。

图像由中途生成，遵循开源许可证。

我们刚刚重新发现了著名的欧拉数，e=2.71828…。结果是，B(x)=e^x。这对我来说并不立即显而易见，但我们可以通过回到方程式 (1) 来看到这一点。

在第二个方程式中我们有e，但在第一个方程式中没有。括号内的x/n 项有点阻碍了这个。为了清理它，让我们通过定义变量来更改：

这将得到方程式 (1)：

注意，像我们上面那样将x取出极限是允许的，适用于连续函数。

指数分布

所以这就是复利和数字e的动机。这一切与点过程和指数分布有什么关系？指数分布在连续时间内工作，并模拟直到某些事件（如车祸）的时间。

理解这一点的最好方法是考虑抛硬币的极限情况。

指数分布的声誉来自于它的无记忆性。事实上，它是唯一的连续无记忆分布。如果你在等待一个到达时间服从指数分布的公交车，那么已经等待了多久并不重要。

你需要等待的额外时间的分布是一样的，无论你刚到达还是已经等了十个小时。这一特性使得指数分布非常容易处理。

当我们将事物离散化时，更容易理解这一属性。与其在连续时间中等待公交车到来，不如想象每分钟抛一次硬币，等待看到正面。我们在看到第一个正面之前看到的反面次数是一个离散随机变量，因为它只能取非负整数值（与公交车到达时间不同，公交车到达时间可以是任何实数，如 3.4 分钟）。这种离散分布称为几何分布。

一个现实世界中的例子是赌场中的老虎机，赌徒不断玩直到赢得大奖。每次旋转都是独立的，这意味着几何分布也是无记忆的。当人们认为一台机器很久没有中奖，所以现在“该中奖了”，他们并没有考虑到过程的无记忆性质，从而陷入了“赌徒谬论”。

我们可以通过抛硬币来模拟机器的每个旋转。硬币正面朝上的概率是p。我们开始抛这枚硬币。我们没有看到正面的概率是多少呢？这意味着我们已经看到k次连续的反面。这个概率是：

现在，我们要转到连续时间。因此，我们将连续的时间线划分为离散部分。硬币抛掷发生在每一个这些离散事件上。每个单位时间t被划分为一个大的离散部分数d。

图：每个时间单位被划分为 d 个间隔，每个间隔的末尾抛一次硬币。

现在我们用T表示有趣事件（硬币正面朝上）发生的时间。为了得到T的分布，我们再次关注其生存函数，即T大于某个数t的概率。我们知道，在这个时间之前必须发生总共[t/d]次抛掷（其中[.]表示取整函数）。例如，如果t=10且每个时间单位被分成3部分，那么到那时已经发生了*[10/3] = [3.33] = 3次抛掷。为了使这成为一个真正的连续时间过程，我们需要使d足够小，以至于消失。但随着我们使d变小，我们会发现抛掷次数增加。因此，我们的p也必须变小以进行补偿（否则，事件将变得非常频繁，以至于任何极小的时间间隔都会有许多事件发生）。因此，p和d变量必须同时趋近于0*。利用上述离散情况的方程，得到时间t之前发生的抛掷次数以及我们的p和d变量必须趋近于0，我们得到T的生存函数：

第二个极限就变成 1。

方程（3）

只有当p和d以线性关系同时减小到0时，这个极限才有趣。因为它们同时趋近于零，直线的截距必须为0。假设直线是：

这是大多数人遇到麻烦的步骤。为什么突然出现这个方程？它从哪里来？如果我们问方程（3）中的双重极限等于什么，答案将是“这要看情况”。取决于p和d之间的关系。首先，我们知道p和d正在一起趋近于0。因此，它们之间的关系必须经过*(0,0)。接下来，我们需要选择它们之间的函数形式。我们可以选择任何形式。但如果我们选择除线性关系以外的任何形式，我们会得到一个平凡的答案（例如任何t下的0或1*），而得不到有趣的连续分布。

通过上述线性关系，方程（3）变为：

方程（3）

只有当p和d以线性关系同时减小到0时，这个极限才有趣。因为它们同时趋近于零，直线的截距必须为0。假设直线是：

上述方程（3）变为：

我们需要另一个替换，以使其与方程（1）对齐：

这是指数分布的生存函数。我们从几何分布出发，取极限并推导出指数分布，同时使用了复利的结果。

理解列表推导式以编写更简洁、更快速的 Python 代码

原文：towardsdatascience.com/comprehending-comprehensions-to-write-cleaner-faster-python-d18908b42c84?source=collection_archive---------19-----------------------#2023-03-20

Charles Mendelson

·

关注 发布于 Towards Data Science ·5 分钟阅读·2023 年 3 月 20 日

–

图片由 Jonathan Cooper 提供，来源于 Unsplash

使用函数和海象运算符编写更好的推导式

TLDR

1. 将转换逻辑封装在函数中。

2. 在列表推导式中调用该函数。

3. 使用海象运算符在列表推导式命名空间内赋值

问题：

列表和字典推导是中级到高级 Python 的基石。不幸的是，它们也很容易被误解为for循环的一对一替代品，当你想要生成一个可迭代对象时。

推导式的语法使得将for循环的内容直接转移到推导式中时，编写清晰的代码变得具有挑战性。

入门：

让我们从一个复杂度适中的问题开始。

对于由随机名称生成器生成的岛屿列表：

1. 移除空白

2. 替换单词“the”

3. 将新单词的首字母变为小写

一个初学者使用for循环可能会这样解决问题：

island_names = ['The Neglected Holm', 'The Barnacle Key', 'The Dancing Enclave',
                'The Peaceful Island', 'Allerport Reef', 'Cresstead Archipelago',
                'Petromeny Islet', 'Esterisle Peninsula', 'Traygami Cay',
                'Savaside Peninsula']

new_island_names_loop = []

for island in island_names:
    new_island = island.replace(' ', '')
    new_island = new_island.replace('The', '')
    # turns first letter into lowercase
    new_island = new_island[0].lower() + new_island[1:]
    new_island_names_loop.append(new_island)

这段代码有效，并且逻辑清晰，但性能和 Pythonic 风格上并不是特别出色。

难以理解的推导式：

仅用列表推导来替换它会产生这段难以理解的乱码：

island_names = ['The Neglected Holm', 'The Barnacle Key', 'The Dancing Enclave',
                'The Peaceful Island', 'Allerport Reef', 'Cresstead Archipelago',
                'Petromeny Islet', 'Esterisle Peninsula', 'Traygami Cay',
                'Savaside Peninsula']

# This is incomprehensible and no one sane would approach the problem this way

new_island_incomprehensible_comprehenson = [
    island.replace(' ', '').replace('The', '')[0].lower() +
    island.replace(' ', '').replace('The', '')[1:]
    for island in island_names]

要了解如何修复这个问题，我们必须理解列表推导的起源。Python 从 Haskell 中借鉴了列表和字典推导，Haskell 是一种非常有主张的函数式语言。要有效地使用它们，你需要以函数式思维来考虑。将可迭代对象中的每一项传递给一个函数来执行转换。

函数式推导

在这段代码中，我将转换逻辑封装在一个函数中，并在推导式中调用该函数，将逻辑与实现分开：

island_names = ['The Neglected Holm', 'The Barnacle Key', 'The Dancing Enclave',
                'The Peaceful Island', 'Allerport Reef', 'Cresstead Archipelago',
                'Petromeny Islet', 'Esterisle Peninsula', 'Traygami Cay',
                'Savaside Peninsula']

# Understanding that list comprehensions come from functional languages
# helped me realize you're supposed to define a function for complex # transformations
# This also helps improve code quality by putting the transformation logic inside
# a function

def prep_island_name(island):
    new_island = island.replace(' ', '')
    new_island = new_island.replace('The', '')
    return new_island[0].lower() + new_island[1:]

new_island_comprehension = [prep_island_name(island) for island in island_names]

在这里，我们保留了原始循环的清晰度，但能够利用推导式的性能提升。

使用for循环也是一种好的实践。将逻辑与实现分开使代码更具模块化。

如果我们需要在推导式中使用相同的函数怎么办？

如果我们的问题需要保留多个转换步骤，并多次使用相同的函数怎么办？

让我们稍微修改一下提示。

现在我们想对我们的岛屿列表进行操作，并：

1. 移除空白

2. 替换单词“the”

3. 将清理后的字符串作为字典中的键

4. 将相同的字符串（首字母小写）作为值附加。

初学者可能会像以前一样使用for循环解决问题：

island_names = ['The Neglected Holm', 'The Barnacle Key', 'The Dancing Enclave', 'The Peaceful Island',
                'Allerport Reef', 'Cresstead Archipelago', 'Petromeny Islet', 'Esterisle Peninsula', 'Traygami Cay',
                'Savaside Peninsula']

new_island_dict = {}
for island in island_names:
    new_island = island.replace(' ', '')
    new_island = new_island.replace('The', '')
    # turns first letter into lowercase
    new_island_dict[new_island] = new_island[0].lower() + new_island[1:]

与我们之前的循环类似，这很容易理解，逻辑清晰，但繁琐。

使用函数式推导

我们将对之前的函数做几个改动。我们将把逻辑分成两个函数，一个用于替换不需要的字符，另一个用于将首字母转换为小写：

island_names = ['The Neglected Holm', 'The Barnacle Key', 'The Dancing Enclave',
                'The Peaceful Island', 'Allerport Reef', 'Cresstead Archipelago',
                'Petromeny Islet', 'Esterisle Peninsula', 'Traygami Cay',
                'Savaside Peninsula']

def prep_island_name(island):
    new_island = island.replace(' ', '')
    new_island = new_island.replace('The', '')
    return new_island

def turn_first_letter_lowercase(string):
    return string[0].lower() + string[1:]

# One challenge of using a dictionary comprehension is that we can end up making
# the same function call multiple times.
new_island_comprehension = {prep_island_name(island): 
turn_first_letter_lowercase(prep_island_name(island))
                            for island in island_names}

这个逻辑是合理的，但我们需要调用prep_island_name()两次，这样显得繁琐且效率低下。

海象操作符来救援

海象操作符:=（也称为赋值操作符），在 Python 3.8 中引入，允许你在传统上不能使用=的命名空间中定义变量名。

现在，除了两次调用 prep_island_name()函数外，我们可以使用赋值运算符在列表推导内部定义一个变量，并将该变量传递给 turn_first_letter_lowercase()。请特别注意我们列表推导中包围关键定义的括号：

island_names = ['The Neglected Holm', 'The Barnacle Key', 'The Dancing Enclave',
                'The Peaceful Island', 'Allerport Reef', 'Cresstead Archipelago',
                'Petromeny Islet', 'Esterisle Peninsula', 'Traygami Cay',
                'Savaside Peninsula']

def prep_island_name(island):
    new_island = island.replace(' ', '')
    new_island = new_island.replace('The', '')
    return new_island

def turn_first_letter_lowercase(string):
    return string[0].lower() + string[1:]

# by using a walrus operator we can call the function once and then assign it a
# variable name in the comprehension name space, avoiding redundant function
# calls
new_island_comprehension_walrus = {
    (new_island := prep_island_name(island)): turn_first_letter_lowercase(new_island)
    for island in island_names}

结论

我花了令人尴尬的长时间才弄明白如何有效地使用列表推导，因为我一直把它们认为是 for 循环。功能性编程帮助我打破了对 for 循环的依赖，并编写了更简洁、更快速的代码。

学习海象运算符去除了 for 循环的最后一个优势——轻松的变量赋值。

关于

查尔斯·门德尔森是一位总部位于西雅图的数据工程师，同时也是华盛顿大学继续教育学院的教学助理，在那里他教授 Python 证书课程。

他即将从哈佛延伸学院获得心理学硕士学位。如果你想与他联系，最好的方式是通过LinkedIn，他曾被 Databand.ai 评选为 2022 年 25 位数据工程影响者之一。

最初发表于 https://charlesmendelson.com 2023 年 3 月 20 日。

Python 中的并发与并行综合指南

原文：towardsdatascience.com/comprehensive-guide-to-concurrency-and-parallelism-in-python-415ee5fbec1a

使用 multiprocessing、threading 和 asyncio

Oliver S

·发表于Towards Data Science ·8 分钟阅读·2023 年 4 月 14 日

–

在这篇文章中，我们将详细介绍 Python 中的并发和并行。我们将介绍这些术语，并展示如何使用multiprocessing、threading和asyncio在 Python 中应用它们。我们将学习何时使用多个进程，何时使用多个线程，并为每种情况提供实际示例。

图片来源于Tomas Sobek的Unsplash

概述

首先，从计算机科学的一般角度介绍标题中的两个术语：并发是指“同时”访问相同资源，例如 CPU 核心或磁盘——而并行则描述多个任务访问不同的资源，例如不同的 CPU 核心。

在 Python 中，通过multiprocessing包实现了并行——它允许创建多个独立的进程。并发可以通过threading包实现，允许创建不同的线程——或者通过asyncio，它遵循稍微不同的哲学。

它们的异同是什么？一个进程是由底层操作系统管理的独立进程。一个进程可以启动多个线程——线程可以被视为进程的一个子例程。默认情况下，进程是独立的实体，不共享内存或类似资源。它们的创建会产生开销，同时数据共享和传递也不是简单的，需要通过进程间通信（IPC）来管理。相比之下，线程是轻量级的，数据共享很容易，因为它们属于同一个进程和内存空间。

Python 与其他几种语言的不同之处在于，它使用了一个全局解释器锁 (GIL)来管理对 Python 解释器的并发访问。这个锁在多线程应用程序中变得有效，确保一次只能有一个线程可以运行 Python 代码。因此，每个多线程 Python 应用程序实际上都是单核的！

有了这些，我们现在可以定义使用场景和建议：何时使用多处理，何时使用多线程：如果您的应用程序是 CPU 绑定的，意味着对 CPU 的访问是主要瓶颈，则使用多处理。只有这样，您的应用程序才能有效地同时使用多个核心，并加速可以在多个核心上并行的代码。缺点是创建进程的开销增加，以及如果数据需要共享则增加的复杂性。然而，除了 CPU 之外，还可能存在其他瓶颈，例如 I/O 绑定的应用程序。这些应用程序的特点是输入/输出操作的长时间等待，例如等待用户输入或等待 web 请求返回。如果是这种情况，多线程可能是更好的选择，因为我们可以避免多处理带来的开销，并原生共享数据。

接下来，我们将介绍这些概念，从多处理开始。然后，我们将使用多个线程实现相同的 CPU 绑定示例应用程序，并显示它确实一次只能限制一个核心。接着，我们给出一个使用 I/O 绑定应用程序的示例，并使用threading和asyncio来实现。

CPU 绑定应用程序

为了在 Python 中实现并行，我们使用 multiprocessing 包。通过使用这个包，我们可以原生定义进程，通过 Process 类，然后简单地启动和停止它们。以下示例启动四个进程，每个进程都计算到 100000000。这意味着，该应用程序是 CPU 绑定的——CPU（核心）增量计数的速度越快，完成的速度就越快：

from multiprocessing import Process
import time

MAX_COUNT = 100000000
NUM_PROCESSES = 4

def count(max_count: int) -> int:
    counter = 0
    for _ in range(max_count):
        counter += 1
    print("Finished")

if __name__ == "__main__":
    start_time = time.time()

    processes = [Process(target=count, args=(MAX_COUNT,)) for _ in range(NUM_PROCESSES)]

    for process in processes:
        process.start()

    for process in processes:
        process.join()

    print(f"Time elapsed: {time.time() - start_time}")

在我的笔记本电脑上，上面的程序执行大约需要 3 秒。

请注意，__main__ 的检查在这种情况下很重要，因为新进程将基于相同的代码生成。如果没有检查，我们会遇到错误，因为这会触发无限循环。

多线程

现在，让我们看看如何使用多线程来解决相同的问题。为此，我们可以简单地用 threading.Thread 替换 multiprocessing.Process，其他方面的代码基本上是类似的：

import threading
import time

MAX_COUNT = 100000000
NUM_PROCESSES = 4

def count(max_count):
    counter = 0
    for _ in range(max_count):
        counter += 1
    print("Finished")

if __name__ == "__main__":
    start_time = time.time()

    threads = [
        threading.Thread(target=count, args=(MAX_COUNT,)) for _ in range(NUM_PROCESSES)
    ]

    for thread in threads:
        thread.start()

    for thread in threads:
        thread.join()

    print(f"Time elapsed: {time.time() - start_time}")

执行时，这段代码大约需要 10 秒才能完成——比使用不同进程时慢 3 倍。我们可以通过经验观察到 GIL 的作用，确实将并行代码执行限制在一个核心上，正如前面所述。

但现在，让我们稍微深入了解一下多处理，特别是使用一些更高层次的抽象来简化启动多个进程。

multiprocessing.Pool

一种这样的抽象是 multiprocessing.Pool。这是一个方便的函数，用于生成一个工人/进程池，它会自动将给定的工作分配给彼此并执行：

from multiprocessing import Pool

MAX_COUNT = 100000000
NUM_PROCESSES = 4

def count(max_count: int) -> int:
    counter = 0
    for _ in range(max_count):
        counter += 1
    print("Finished")
    return counter

if __name__ == "__main__":
    results = Pool(NUM_PROCESSES).map(count, [MAX_COUNT] * 5)
    print(results)

如我们所见，我们实例化一个 Pool，并设置要使用的工作者数量。建议将此数量设置为您拥有的 CPU 核心数（os.cpu_count）。然后我们使用 map 函数，将我们希望每个工作者运行的函数作为第一个参数传递，将输入数据作为第二个参数传递。在我们的案例中，这只是参数 max_count。由于我们传递了一个长度为 5 的数组，count 将运行 5 次。将池工作者数量设置为 4，结果是 2 个“周期”：在第一个周期中，工作者 0–3 处理前 4 个参数/数据集，在第二个周期中，完成得最早的工作者处理最后一个数据集。在这个示例中，count 还返回一个值，以展示 map 如何处理返回值。

Pool.imap

为了结束这一部分，让我们看看一个非常类似于 map 的函数，即 imap。map 和 imap 之间的区别在于任务处理的方式：map 将所有任务转换为一个列表，然后一次性将它们传递给工作者，将任务分解成块。另一方面，imap 逐个传递任务。这可能更节省内存，但也可能较慢。另一个区别是 imap 在每个任务完成后立即返回，而 map 会阻塞直到所有任务完成。我们可以这样使用 imap：

results = Pool(NUM_PROCESSES).imap(count, [MAX_COUNT] * 5)
for result in results:
    print(result)

I/O 绑定的应用程序

这总结了我们对并行性的介绍，它对于 CPU 绑定的应用程序非常有用并且推荐使用。现在，让我们谈谈那些不是 CPU 绑定的应用程序，特别是 I/O 绑定的应用程序。为此，我们使用并发。特别是，我们将通过一个示例来介绍这个概念，使用多线程，因为对于这种情况，使用多进程没有意义，因为：

• 这个应用程序不是 CPU 绑定的，因此多进程带来的额外开销不值得。

• 线程之间的通信会引入额外的开销和复杂性，当与多进程一起使用时尤为如此。

我们选择了以下示例来表示一个 I/O 绑定的应用程序：有一个发布者线程生成数据，在我们的示例中，这些数据来自用户。然后将有 N 个订阅者线程，它们等待数据中的特定条件，一旦这个条件激活，就开始执行一些涉及大量空闲操作的任务，即 CPU 不是瓶颈。

我们可以看到，对于这种情况，CPU 不是限制因素：大量时间将花费在等待输入上，并且随后的操作也不是 CPU 密集型的。然而，我们仍然希望/需要某种形式的并发——我们有不同的线程独立运行（生产者与订阅者），并且还希望执行后续任务时能够“并行”进行，而不是顺序执行。

因此，我们使用 threading 包，示例代码如下：

import threading
import time

NUM_CONSUMERS = 2

condition_satisfied = False

lock = threading.Lock()

def producer() -> None:
    global condition_satisfied

    while True:
        user_input = input("Enter a comamnd:")
        if user_input == "Start":
            # Signal the producers to start
            lock.acquire()
            condition_satisfied = True
            lock.release()
            break
        else:
            print(f"Unknown command {user_input}")
        lock.release()
        time.sleep(1)

def consumer(consumer_idx: int) -> None:
    while True:
        lock.acquire()
        condition_satisfied_read = condition_satisfied
        lock.release()
        if condition_satisfied_read:
            for i in range(10):
                print(f"{i} from consumer {consumer_idx}")
                time.sleep(1)
            break
        time.sleep(1)

if __name__ == "__main__":
    threads = [threading.Thread(target=producer)] + [
        threading.Thread(target=consumer, args=(consumer_idx,))
        for consumer_idx in range(NUM_CONSUMERS)
    ]

    for thread in threads:
        thread.start()

    for thread in threads:
        thread.join()

生产者线程不断查询用户输入，一旦用户输入“开始”，一个标志会被切换，通知消费者开始。在消费者中，我们对此做出反应，并在收到信号后在 10 秒内计数到 10。

asyncio

为了总结这篇文章，我们想快速介绍一下 [asyncio](https://medium.com/towards-data-science/introduction-to-asyncio-57a5a1290ce0) 以及如何在这种情况下使用它。对于完整的介绍，我建议参考链接中的文章。在这里，我们简单说一下 asyncio 可以被视为一种稍微不同的多线程代码编写/结构化范式——许多开发者在多种场景中更喜欢这种方式，主要是因为它的简洁性。asyncio 适用于 I/O 密集型应用，它采用了协程的概念：当一个协程在等待某个条件（如用户输入、网络请求完成等）时，它可以“让出”控制权，让另一个协程运行。

使用这个原则，上面的示例可以这样编写：

import asyncio

NUM_CONSUMERS = 2

condition_satisfied = False
should_terminate = False

async def producer():
    global condition_satisfied

    while True:
        user_input = input("Enter a comamnd:")
        if user_input == "Start":
            # Signal the producers to start
            condition_satisfied = True
            break
        else:
            print(f"Unknown command {user_input}")
        await asyncio.sleep(1)

async def consumer(consumer_idx):
    global condition_satisfied

    while True:
        if condition_satisfied:
            for i in range(10):
                print(f"{i} from consumer {consumer_idx}")
                await asyncio.sleep(1)
            break
        await asyncio.sleep(1)

async def main():
    await asyncio.gather(producer(), consumer(0), consumer(1))

if __name__ == "__main__":
    asyncio.run(main())

请注意，在这种情况下不需要使用锁，因为我们决定何时让出执行并切换到不同的线程：我们仅在“等待”期间进行切换，这时没有共享变量被写入或访问。

结论

在这篇文章中，我们介绍了并行性和并发性概念，并描述了它们如何转换为 Python 中的多进程和多线程。

由于 GIL，Python 中的多线程应用实际上是单核的——因此，这种范式不适用于 CPU 密集型应用。对于这些应用，我们展示了如何管理多进程——然后通过将相同的示例重写为使用多个线程，经验性地证明了由于 GIL 导致的性能下降。

尽管如此，仍然有许多场景中多线程也是有益和必要的，特别是对于 I/O 密集型应用。对于这些应用，我们介绍了 Python 中的多线程，并通过将示例转换为 asyncio 这一稍微不同的多线程应用编写范式来结束。

这结束了本教程，希望对你有所帮助。下次见！

排名评估指标的综合指南

原文：towardsdatascience.com/comprehensive-guide-to-ranking-evaluation-metrics-7d10382c1025?source=collection_archive---------0-----------------------#2023-07-29

探索丰富的指标选择，并找到适合你问题的最佳指标

Vyacheslav Efimov

·

关注 发表在 Towards Data Science · 13 分钟阅读 · 2023 年 7 月 29 日

–

介绍

排名是机器学习中的一个问题，其目标是以最合适的方式对文档列表进行排序，以便最相关的文档出现在顶部。排名出现在数据科学的多个领域，从推荐系统开始，在这些系统中，算法建议一组购买项目，直到 NLP 搜索引擎，通过给定的查询，系统尝试返回最相关的搜索结果。

自然产生的问题是如何评估排名算法的质量。与经典的机器学习一样，并不存在适用于任何类型任务的单一通用指标。为什么？原因很简单，因为每种指标都有其适用范围，这取决于特定问题的性质和数据特征。

这就是为什么了解所有主要指标对于成功解决任何机器学习问题至关重要。正是我们将在本文中做到的。

然而，在继续之前，让我们了解为什么某些流行指标通常不应被用于排名评估。通过考虑这些信息，将更容易理解其他更复杂指标存在的必要性。

注意。本文及使用的公式基于Ilya Markov 的离线评估演示。

指标

我们将在本文中讨论几种信息检索指标：

不同类型的指标

未排序的指标

设想一个推荐系统预测电影评分并向用户展示最相关的电影。评分通常表示一个正的实数。乍看之下，像MSE（RMSE，MAE等）这样的回归指标似乎是评估系统在保留数据集上的质量的合理选择。

MSE考虑了所有预测的电影，并测量真实标签与预测标签之间的平均平方误差。然而，最终用户通常只对显示在网站首页的前几项结果感兴趣。这表明，他们对搜索结果末尾的低评分电影并不真正感兴趣，这些电影在标准回归指标下也被同等估计。

下面的简单示例演示了一对搜索结果，并测量每个结果中的MSE值。

对两个查询的错误估计表明 MSE 是排名的糟糕指标。绿色文档是相关的，而红色文档是不相关的。文档列表按预测相关性排序（从左到右）。

尽管第二个搜索结果的MSE较低，但用户不会对这样的推荐感到满意。首先查看非相关项时，用户必须向下滚动才能找到第一个相关项。因此，从用户体验的角度来看，第一个搜索结果要好得多：用户对顶部的项感到满意，并继续使用，而不关心其他项。

相同的逻辑也适用于分类指标（精确度，召回率），这些指标也会考虑所有项目。

精确度和召回率公式

所有描述的指标有什么共同点？它们都平等对待所有项，并且不考虑高相关和低相关结果之间的差异。这就是为什么它们被称为未排序的原因。

通过上述两个相似的问题示例，我们在设计排名指标时需要关注的方面似乎更加明确：

一个排名指标应该对更相关的结果给予更高的权重，同时降低或忽略不太相关的结果。

排名指标

Kendall Tau 距离

Kendall Tau 距离 基于排名逆序对的数量。

逆序对是指文档 (i, j) 的一对，其中文档 i 的相关性大于文档 j，但在搜索结果中出现在 j 之后。

Kendall Tau 距离计算排名中的所有逆序对。逆序对越少，搜索结果越好。尽管该指标看起来合乎逻辑，但它仍然有一个缺点，这在下面的示例中有所展示。

尽管逆序对的数量较少，但从用户的角度来看，第二个排名仍然更差

似乎第二个搜索结果只有 8 个逆序对，比第一个结果的 9 个要好。类似于上面的MSE示例，用户只对第一个相关结果感兴趣。在第二种情况下经过几个不相关的搜索结果，用户体验会比第一种情况更差。

Precision@k & Recall@k

与通常的精确度和召回率不同，我们可以考虑仅关注某个数量的前* k 个推荐项。这样，指标就不关心低排名的结果。根据所选的 k 值，相应的指标被称为precision@k*（“k 时的精确度”）和recall@k（“k 时的召回率”）。它们的公式如下所示。

precision@k 和 recall@k 的公式

想象一下，前* k *个结果显示给用户，其中每个结果可能相关也可能不相关。precision@k 衡量前 k 个结果中相关结果的百分比。同时，recall@k 评估前 k 个相关结果占整个数据集中相关项总数的比例。

为了更好地理解这些指标的计算过程，让我们参考下面的示例。

precision@k 和 recall@k 计算示例。绿色文档表示相关项，而红色文档表示不相关项。

系统中有 7 个文档（命名为 A 到 G）。根据其预测，算法从中选择 k = 5 个文档给用户。如我们所见，前 k = 5 中有 3 个相关文档 (A, C, G)，这导致 precision@5 等于 3 / 5。同时，recall@5 考虑了整个数据集中相关项目的数量：其中有 4 个 (A, C, F 和 G)，使得 recall@5 = 3 / 4。

recall@k 随着 k 的增长总是增加，这使得该指标在某些情况下不完全客观。在边缘情况下，当系统中的所有项目都显示给用户时，recall@k 的值等于 100%。precision@k 不具备与 recall@k 相同的单调性，因为它衡量的是与前 k 个结果相关的排名质量，而不是整个系统中相关项目的数量。客观性是 precision@k 通常被实际应用中优于 recall@k 的原因之一。

AP@k (平均精度) & MAP@k (均值平均精度)

普通的 precision@k 问题在于它没有考虑相关项在检索文档中的出现顺序。例如，如果检索到 10 个文档，其中 2 个是相关的，precision@10 将始终相同，尽管这 2 个文档在 10 个文档中的位置不同。例如，如果相关项位于位置 (1, 2) 或 (9, 10)，该指标会区分这两种情况，结果 precision@10 为 0.2。

然而，在实际应用中，系统应该对排名靠前的相关文档给予更高的权重，而不是排名靠后的。这一问题通过另一种指标 平均精度 (AP) 得以解决。与普通 precision 相同，AP 的取值范围在 0 和 1 之间。

平均精度公式

AP@k 计算从 1 到 k 所有值的 precision@i 的平均值，其中第 i 个文档是相关的。

对两个查询计算的平均精度

在上图中，我们可以看到相同的 7 个文档。对查询 Q₁ 的响应产生了 k = 5 个检索文档，其中 3 个相关文档位于索引 (1, 3, 4)。对每个位置 i，计算 precision@i：

• precision@1 = 1 / 1

• precision@3 = 2 / 3

• precision@4 = 3 / 4

所有其他不匹配的索引 i 被忽略。AP@5 的最终值是以上精度的平均值：

• AP@5 = (precision@1 + precision@3 + precision@4) / 3 = 0.81

为了比较，让我们看看对另一个查询 Q₂ 的响应，该查询在前 k 个结果中也包含 3 个相关文档。然而，这次，2 个不相关的文档位于比前一个情况更高的位置（在位置 (1, 3))，这导致 AP@5 较低，为 0.53。

有时需要评估算法在多个查询上的质量，而不是单个查询。为此，使用了均值平均精度（MAP）。它简单地取多个查询 Q 中 AP 的平均值：

平均精度公式

以下示例展示了如何计算针对 3 个不同查询的MAP：

计算三个查询的 AP 和 MAP

RR（倒数排名）和 MRR（均值倒数排名）

有时用户只对第一个相关结果感兴趣。倒数排名是一种度量，返回一个介于 0 和 1 之间的数字，指示第一个相关结果距离顶部的远近：如果文档位于位置 k，则RR的值为 1 / k。

与AP和MAP类似，均值倒数排名 (MRR) 衡量多个查询中的平均RR。

RR 和 MRR 公式

以下示例展示了如何计算 3 个查询的RR和MRR：

三个查询计算出的 RR 和 MRR

用户导向的度量

尽管排名度量考虑了项的排名位置，因此相对于未排名的度量更为优选，但它们仍然有一个显著的缺点：没有考虑用户行为的信息。

以用户为导向的方法对用户行为做出特定假设，并基于这些假设生成更适合排序问题的度量。

DCG（折扣累积增益）和 nDCG（归一化折扣累积增益）

DCG 度量的使用基于以下假设：

高度相关的文档在搜索引擎结果列表中越早出现（排名越高）就越有用 — 维基百科

这一假设自然表示用户如何评估较高的搜索结果，与较低的结果相比。

在DCG中，每个文档都被分配一个增益值，表示特定文档的相关性。对于每个项，给定真实的相关性 Rᵢ（实际值），存在几种定义增益的方法。其中一种最流行的方法是：

DCG 中可能的增益公式

基本上，指数对相关项施加了强烈的强调。例如，如果给电影的评分在 0 到 5 之间，那么每部电影的相应评分将大约具有双倍的重要性，相比于评分减少 1 的电影：

相关性的增益函数

除此之外，根据排名位置，每个项会接收一个折扣值：项的排名位置越高，对应的折扣就越高。折扣作为一种惩罚，通过按比例减少项的增益。在实际应用中，折扣通常选择作为排名指数的对数函数：

DCG 中的折扣公式

排名位置的折扣函数

最后，DCG@k 定义为前 k 个检索项的增益与折扣之和：

一般的 DCG 公式

将 gainᵢ 和 discountᵢ 替换为上述公式后，表达式变为以下形式：

DCG 公式

为了使 DCG 指标更具可解释性，它通常通过完美排名时的 DCGₘₐₓ 最大可能值进行归一化，其中所有项按照其相关性被正确排序。结果的指标称为 nDCG，其值介于 0 和 1 之间。

nDCG 公式

下图展示了 5 个文档的 DCG 和 nDCG 计算示例。

对一组检索文档计算的 DCG 和 nDCG

RBP（排名偏倚精度）

在 RBP 工作流中，用户并不打算检查所有可能的项。相反，他或她以概率 p 从一个文档顺序地移动到另一个文档，并以相反的概率 1 — p 在当前文档处终止搜索程序。每次终止决策都是独立的，不依赖于搜索的深度。根据所做的研究，这种用户行为在许多实验中被观察到。根据 Rank-Biased Precision for Measurement of Retrieval Effectiveness 中的信息，工作流可以在下图中说明。

参数 p 被称为 持久性。

RBP 模型工作流

在这一范式中，用户总是查看第 1 个文档，然后以概率 p 查看第 2 个文档，以概率 p² 查看第 3 个文档，依此类推。最终，查看文档 i 的概率变为：

用户只有在文档 i 已经被查看并且搜索程序立即以概率 1 — p 终止时才会检查文档 i。

之后，可以估算期望的已查看文档数量。由于 0 ≤ p ≤ 1，下列级数是收敛的，表达式可以转化为以下格式：

同样，给定每个文档的相关性 Rᵢ，我们可以找到预期的文档相关性。更高的预期相关性值表示用户对他或她决定查看的文档会更满意。

最后，RPB 计算为期望文档相关性（效用）与期望检查文档数量的比率：

RPB 公式确保其值在 0 和 1 之间。通常，相关性分数是二值型（如果文档相关则为 1，否则为 0），但也可以取 0 到 1 之间的实际值。

应根据用户在系统中的持久性来选择适当的 p 值。较小的 p 值（小于 0.5）更强调排名靠前的文档。较大的 p 值则减少了对前几个位置的权重，并将其分配到较低的位置。有时可能很难找到一个好的持久性 p 值，因此最好进行几次实验并选择效果最佳的 p。

ERR（期望倒数排名）

正如名称所示，这一指标测量了许多查询的平均倒数排名。

该模型类似于 RPB，但有一点不同：如果当前项目对用户是相关的（Rᵢ），则搜索过程结束。否则，如果项目不相关（1 — Rᵢ），则用户以概率 p 决定是否继续搜索过程。如果是，则搜索继续到下一个项目。否则，用户结束搜索过程。

ERR 模型工作流程

根据 Ilya Markov 的离线评估演示，我们来找出 ERR 计算的公式。

首先，让我们计算用户查看文档 i 的概率。基本上，这意味着所有 i — 1 个之前的文档都不相关，并且在每次迭代中，用户以概率 p 继续查看下一个项目：

如果用户停在文档 i 上，这意味着该文档已经被查看，并且以概率 Rᵢ，用户决定终止搜索过程。与此事件相对应的概率实际上与倒数排名相同，即 1 / i。

从现在开始，简单地使用期望值公式，可以估计期望的倒数排名：

参数 p 通常选择接近 1。

与 RBP 的情况一样，Rᵢ 的值可以是二值的，也可以是 0 到 1 之间的实际值。下面的图示例演示了 ERR 计算的一个例子，涉及一组 6 个文档。

ERR 计算。左侧和右侧分别显示了最佳和最差的可能排名。为简化起见，参数 p 选择为 1。

左侧，所有检索到的文档按相关性降序排列，得出最佳的ERR。与右侧情况相反，文档按相关性升序排列，导致最差的ERR。

ERR 公式假设所有相关性分数的范围是 0 到 1。如果初始相关性分数超出这个范围，则需要进行归一化。最常见的做法之一是通过指数归一化：

结论

我们讨论了信息检索中用于质量评估的所有主要指标。用户导向的指标更常使用，因为它们反映了真实用户行为。此外，nDCG、BPR和ERR指标相对于我们迄今为止看到的其他指标具有优势：它们能够处理多个相关性级别，使其更具通用性，而像AP、MAP或MRR这样的指标仅适用于二元相关性级别。

不幸的是，所有这些描述的指标要么是不连续的，要么是平坦的，使得在问题点的梯度为 0 或甚至未定义。因此，大多数排名算法难以直接优化这些指标。然而，许多研究已经在这一领域展开，许多先进的启发式方法已经出现在最流行的排名算法的背后，以解决这个问题。

资源

• 肯德尔τ距离 | 维基百科

• 排名偏差精度测量检索效果

• 折扣累积增益 | 维基百科

• 排名偏差精度的不确定性

• 信息检索，离线评估 | 伊利亚·马尔科夫

除非另有说明，所有图像均为作者所作。

综合时间序列探索性分析

原文：towardsdatascience.com/comprehensive-time-series-exploratory-analysis-78bf40d16083

深入分析空气质量数据

Erich Silva

·发表于Towards Data Science ·阅读时间 18 分钟·2023 年 11 月 25 日

–

由Jason Blackeye拍摄的照片，刊登在Unsplash

你面临一个按时间戳索引的数据集。你的数据可能涉及存储需求和供应，你的任务是预测一个战略产品的理想补货间隔。或者你需要将历史销售信息转化为团队的关键行动洞察。也许你的数据是财务数据，包括历史利率和一系列股票价格。也许你需要建模市场波动性，并量化投资期限内的货币风险。从社会科学和能源分配到医疗保健和环境研究，例子数不胜数。但这些场景有什么共同点？第一，你有一个时间序列任务。第二，你肯定会从开始时进行简明而全面的探索性分析中受益。

目录

• 本文的目标

• 数据集描述

• 库和依赖

• 开始使用

• 全局视角

• 详细视图

• 缺失值

• 间歇性

• 季节性

• 皮尔逊相关

• 平稳性

• 一阶差分

• 自相关

• 参考文献

本文的目标

但进行探索性时间序列分析意味着什么呢？与其他数据科学问题不同，从时间序列数据中获取见解可能是棘手的，一切都不是简单的。你的数据可能具有重要的潜在趋势和季节性，或适合于其复杂的周期模式中的嵌套预测。区分由于数据生成过程中的失败造成的异常离群值与实际的异常情况（它们包含关键信息）可能是具有挑战性的。处理缺失值也可能不像你预期的那么简单。

本文将概述一个在研究时间序列数据集时对我有效的过程。你将跟随我探索细颗粒物的测量值，即 PM 2.5，它是空气污染和空气质量指数的主要来源之一。我将专注于制定一些最佳实践，特别注意细节，以生成清晰且高度信息化的可视化和统计摘要。

数据集描述

这里研究的数据来自加拿大不列颠哥伦比亚省温哥华市的四个监测站。它们包含从 2016 年 1 月 1 日到 2022 年 7 月 3 日的 PM 2.5（直径为 2.5 微米及以下的细颗粒）的一小时平均测量值，单位为 µg/m3（每立方米微克）。

PM 2.5 主要来自燃烧化石燃料，在城市中，通常源自汽车交通和建筑工地。另一个主要的污染源是森林和草地火灾，它们很容易被风吹走 [1]。

下图显示了我们将要探索的站点的大致位置。

图 1. 温哥华地图及空气监测站。作者在 Google Maps 中定制的地图。

数据集来自 British Columbia Data Catalogue，根据发布者的说法，它尚未经过质量保证 [5]。对于你将在这里看到的版本，我已经处理了一些小问题，例如将负值测量（仅 6 个观察值中的 57k）标记为缺失值，并生成了一个包含我们选择的站点的主 DataFrame。

库和依赖

我们将使用 Python 3.9 以及绘图库 Matplotlib 和 Seaborn 来进行可视化。对于统计测试和数据探索，我们将使用 statsmodels Python 模块和 SciPy 库。所有的数据处理和辅助任务将使用 Pandas 和 Numpy 处理。

这些软件包在流行的 Python 发行版和托管的笔记本中原生提供，如 Anaconda 和 Miniconda、Google Colab 或 Kaggle Notebooks。因此，本文中的每个代码示例都应该在你选择的环境中轻松复现。

开始使用

从导入我们的库开始，我们将调用matplotlib.dates.mdates和datetime模块来帮助我们处理 DateTime 索引。为了生成一致的可视化，我还喜欢先定义我们的图表样式和颜色调色板。那么让我们开始吧。

图 2. Seaborn “mako” 颜色调色板。图片由作者提供。

在读取.csv 文件后，我们将定义时间戳DATE_PS为 NumPy 的datetime64对象，并将其设为我们的 DataFrame 索引。这一步将启用一些 Pandas 时间序列功能，例如下文中用于在数据集中创建日期部分特征的功能。

图 3. 带有日期部分特征的主 DataFrame 切片。图片由作者提供。

全景视图

这里是我们将深入探讨的广阔视图——这是我们花一些时间初步了解数据的地方。

对于这个可视化，我们将使用 Seaborn 的关系图，它将从我们 DataFrame 的聚合长格式中读取。为此，我们将使用 Pandas 的melt和resample方法，并在 24 小时间隔内进行均值聚合。这将把数据粒度从每小时减少到每日平均测量值，从而减少生成图表所需的时间。

图 4. 监测站 PM 2.5 时间序列图。图片由作者提供。

通过在整个时间序列跨度内清晰地了解所有四个车站的数据，我们已经可以开始做一些记录：

• 存在一些主要的异常情况，它们似乎在夏季和初秋期间较为普遍。

• 这些异常似乎是由于大规模事件引起的，因为它们影响了所有四个车站，且大致发生在相同的时间段内。

• 如果仔细观察，我在每个图表中都包含了一张淡灰色的散点图，显示了所有车站的数据。通过这个细微的细节，可以看到例如 2017 年的异常在北温哥华的两个车站产生了重大影响（它们的 PM 2.5 值更高），而 2018 年的情况则相反。这种技术还确保所有四个折线图都在相同的 Y 轴范围内。

从第一个图表中可以汲取一些好的实践：

• Matplotlib 允许高度自定义的轴刻度定位器和格式化器。在这个例子中，我使用了mdates模块中的MonthLocator来创建 X 轴上的小月定位器，这有助于整体可读性。

• 我们的图表标题或副标题中传递了从图表中预期的准确日期范围。之所以用“预期”是因为我们的可视化可能由于绘图周期两端的缺失值而被截断，这种做法可以帮助识别这些问题。这也是一个好的文档实践，记录你的发现。

很好的开始，但让我们稍微缩小视图。在下一节中，我们将查看较短的时间段，但现在保持原有的每小时粒度。

详细视图

从现在开始，我们将开始定义一些可以快速调用以生成定制可视化的函数。你可以将其视为建立一个分析工具集，这将在未来非常有帮助。

第一个函数将帮助我们查看特定时间段内的单独时间序列。我们将从查看 2017 年的北温哥华马洪公园站开始。

图 5. 2017 年北温哥华马洪公园 PM 2.5 图。图像由作者提供。

我们可以看到，除了主要异常外，还有局部较小的波动。在 2017 年年初和 12 月，也有较高波动性（在短时间内方差增加）的时期。

让我们进一步深入，查看异常范围之外的情况，以便在 Y 轴上分析我们在更窄范围内的值。

图 6. 2017 年 4 月 15 日至 7 月 1 日的北温哥华马洪公园 PM 2.5 图。图像由作者提供。

我们可以看到这里有一些缺失值。我们函数的fill=True参数有助于识别这些缺失值，并且是给数据中缺失情况提供视觉重点的好方法。这些最初难以察觉的小中断现在变得非常明显。

另一个你可能注意到的细节是 X 轴的自定义日期格式。对于上面的图，我增强了我们的plot_sequence()函数，添加了自定义的主要和次要定位器及格式化程序。这个新功能现在使我们的图表适应可视化的时间跨度，并相应地格式化 X 轴。下面是包含在函数中的代码片段。

现在我们知道我们的数据集有中断，所以让我们更好地查看一下。

缺失值

对于表格数据问题，在本节中，我们可能会专注于将 MAR（随机缺失）与 MNAR（非随机缺失）区分开来。相反，了解我们数据的性质（传感器时间测量），我们知道数据流中的中断可能不是有意的。因此，在这些情况下，区分孤立的缺失值与连续缺失值，以及完全缺失样本中的缺失值是更重要的。这里的可能性很广泛，如果你想了解更多，我专门撰写了一篇文章：

## 时间序列中的缺口处理

短序列和长序列插补的缺失分析和评估方法

towardsdatascience.com

现在，我们先从查看缺失值热图开始。我们将为此定义一个函数。

图 7. 缺失值热图。图片由作者提供。

热图很有用，因为它们允许我们量化缺失情况并在时间轴上本地化它们。从这里，我们可以标记：

• 我们没有完全缺失的样本（缺失值在一个时间段内同时出现）。这是可以预期的，因为来自监测站的数据流是独立的。

• 在时间轴的早期有很长的缺失值序列，随着时间的推移，数据的可用性似乎有所改善。

在后续部分的一些统计分析中，缺失值会带来问题。因此，我们将使用简单的技术来处理这些缺失值，比如：

• Pandas ffill() 和 bfill() 方法。它们分别用于向前或向后填充最近的可用值。

• 使用 Pandas interpolate() 方法进行线性或样条插值。它利用相邻观察值绘制曲线来填充缺失的区间。

间歇性

根据我们数据的特性，我们不应期望出现负值。如开头所述，我在数据预处理时将其视为缺失值。让我们调用汇总统计来确认这一点。

图 8. 汇总统计。图片由作者提供。

我们看到每个站点的最小测量值为零，这引出了下一个问题。我们的时间序列是否间歇性？

当数据中有大量值恰好为零时，这种现象被称为间歇性。这种行为带来了特定的挑战，必须在模型选择过程中加以考虑。那么，我们的数据中零值出现的频率有多高呢？

图 9. 零值计数。图片由作者提供。

我们可以看到零值的数量可以忽略不计，因此我们没有间歇性序列。这是一个简单但至关重要的检查，特别是如果你的目标是预测的话。对于某些模型来说，预测绝对零可能比较困难，例如，如果你想预测需求，这可能会成为一个问题。你不想计划向你的客户交付三件产品，如果实际上他什么都不期望的话。

季节性

了解时间序列中的周期对规划建模阶段至关重要。如果你决定过度聚合数据，你可能会丢失较小周期的关键信息，或者这可以帮助你确定在更小粒度下进行预测的可行性。

我们将使用一些箱线图来开始查看这些数据。但首先，我们会暂时移除前 5%的百分位数，以便以更好的尺度查看数据。

在下一个函数中，我们将使用一系列箱线图来调查数据中的周期。我们还将把我们的颜色调色板映射到中位数值，以便作为另一个精美的视觉提示。

图 10. PM 2.5 每小时值。图片由作者提供。

这个初步图表返回每小时的测量值。这里我们可以看到：

• 从早上 9 点到下午 2 点，PM 2.5 值 consistently higher。

• 北温哥华以外的站点也在晚上 8 点到 10 点出现峰值。

• 早晨 2 点到 5 点的 PM 2.5 值最低。

现在查看周季节性以及一周内的值差异。

图 11。PM 2.5 每日值。图片由作者提供。

从这里我们可以看到：

• 周末的 PM 2.5 值较低。

• 星期二污染水平有略微上升的趋势。

最后，查看月度趋势。

图 12。PM 2.5 月度值。图片由作者提供。

我们可以观察到：

• 所有年份中 8 月的 PM 2.5 值 consistently higher。

• 南部站点在 6 月和 7 月的 PM 2.5 值较低，而北温哥华的站点在 1 月显示出较低的测量值。

最后，从这些图表中得出的更多良好实践：

• 不要简单使用你的颜色调色板，因为它们可能会误导你得出相同的解释。如果我们只是将pallette=”mako”传递给我们的箱线图，它会映射到 X 轴，而不是我们感兴趣的变量。

• 网格图是低维数据的强大信息容器，可以通过 Seaborn 的relplot()或 Matplotlib 的subplots()快速设置。

• 你可以利用 Seaborn 的boxplot() order参数来重新排序你的 X 轴。我用它重新排列了周几的 X 标签，使之有意义。

更详细的季节性视图可以通过我们时间序列的趋势-季节分解获得。然而，这将留到未来的文章中，我们可以更深入地探讨时间序列相似性和模型选择。如果你对此感兴趣，确保关注我在 Medium 上的未来出版物。

现在，让我们快速查看我们熟知的统计系数，以调查四个站点之间的线性关系。

Pearson 相关

R 程序员可能对以下图表很熟悉。一个相关图是一种简洁且信息丰富的可视化，已在多个 R 库中实现，例如GGally包中的ggpairs()。相关图的上对角线显示了双变量相关性，即我们数据中数值变量之间的 Pearson 相关系数。在下对角线中，我们可以看到带有回归曲线的散点图。最后，在主对角线上，我们展示了每个变量的直方图和密度曲线。

以下代码是使用 Seaborn PairGrid()图表的改编实现，并且是我们分析工具集中的另一个函数。

图 13。所有四个站点的 PM 2.5 相关图。图片由作者提供。

正如预期的那样，我们的站点高度相关，尤其是那些彼此更近的站点，例如北温哥华的两个站点。

重要的是要注意，为了减轻计算时间，我们的数据被聚合成 6 小时的时间段。如果你在这个图上尝试更大的聚合时间段，你会看到相关系数的增加，因为均值聚合往往会平滑数据中的异常值。

如果你已经了解了时间序列分析，你现在可能会考虑其他值得检查的相关性。但首先，我们需要测试我们的时间序列是否平稳。

平稳性

平稳时间序列是指其统计特性随时间不变。换句话说，它具有一个恒定的均值、方差和自相关性，与时间无关[4]。

多个预测模型依赖于时间序列平稳性，因此在这一探索阶段测试其平稳性至关重要。我们的下一个函数将利用 statsmodels 实现的两个常用平稳性测试，即Augmented Dickey-Fuller（“ADF”）测试和Kwiatkowski-Phillips-Schmidt-Shin（“KPSS”）测试。

我将把两个测试的假设列在下面。请注意，它们有相反的原假设，因此我们将创建一个“决策”列，以便于解释它们的结果。你可以在statsmodels 文档中了解更多关于这两个测试的信息。

Augmented Dickey-Fuller (ADF) 测试假设：

• H0: 时间序列样本中存在单位根 (非平稳)

• Ha: 时间序列样本中不存在单位根 (平稳)

Kwiatkowski-Phillips-Schmidt-Shin (KPSS) 测试假设：

• H0: 数据围绕常数是平稳的 (平稳)

• Ha: 时间序列样本中存在单位根 (非平稳)

现在一个恰当的问题是我们应该在什么尺度上检查平稳性。答案将高度依赖于你如何建模数据，而全面的探索性分析的一个目标正是帮助你做出这个决定。

为了说明目的，在接下来的例子中，我们将查看 2016 年 1 月和 2022 年 1 月的温哥华国际机场站数据，看看数据在 2016 年至 2022 年间是否有行为上的变化。

你可能还记得我们在“缺失值”部分提到过，我们可以使用 Pandas 的ffill()、bfill()和interpolate()方法来快速填补序列中的中断。你可以看到，我在我们的函数中定义了一个专用的fillna参数，以便选择这几种方法中的任意一种，快速处理缺失值，因为这两个测试只接受完整样本。

现在回到我们的结果。

图 14. 2016 年 1 月 ADF 和 KPSS 平稳性测试结果。图片由作者提供。

图 15. 2022 年 1 月 ADF 和 KPSS 平稳性测试结果。图片由作者提供。

我们可以看到，2016 年的两个测试都指示非平稳性，但 2022 年的结果则有所不同。statsmodels 文档 清楚地列出了当 ADF 和 KPSS 测试一起进行时的结果解释[6]：

• 案例 1：两个测试均得出系列不是平稳的结论 — 该系列不是平稳的

• 案例 2：两个测试均得出系列是平稳的结论 — 该系列是平稳的

• 案例 3：KPSS 指示平稳性 和 ADF 指示非平稳性 — 该系列是趋势平稳的。需要去除趋势 以使系列严格平稳。去趋势后的系列平稳性进行了检查。

• 案例 4：KPSS 指示非平稳性 和 ADF 指示平稳性 — 该系列是差分平稳的。需要使用差分 使系列平稳。差分系列的平稳性进行了检查。

如果你对所有四个站点进行多个月份的重复操作，你会发现 案例 4 在数据中占主导地位。这将引导我们进入下一个部分，即 使用一阶差分使数据平稳。

一阶差分

作为最常见的转换技术之一，对时间序列应用一阶或二阶差分被广泛使用，以使数据适合只能用于平稳时间序列的统计模型。在这里，我们将查看该技术在 2016 年 1 月的一个示例中的应用。但首先，让我们用我们的 plot_sequence() 函数查看转换前的原始数据。

图 16. 温哥华国际机场 2016 年 1 月的 PM 2.5 图。图片由作者提供。

我们可以看到，期间方差从月初到月底显著变化。均值 PM 2.5 也似乎从较高值变为较低且更稳定的值。这些是确认系列非平稳性的特征之一。

再次，Pandas 提供了一个非常方便的方法来对数据进行差分。我们将调用 .diff() 方法来创建数据的第一阶差分版本。然后再次绘制相同的时间段。

图 17. 温哥华国际机场 2016 年 1 月的 PM 2.5 差分图。图片由作者提供。

除了仍然存在的波动方差外，数据现在在均值附近明显更加稳定。我们可以再次调用我们的 stationarity_test() 函数来检查差分数据的平稳性。

图 18. 差分数据的 ADF 和 KPSS 平稳性测试结果。图片由作者提供。

这就是结果。我们可以对我们全面的探索性时间序列分析进行另一项检查，因为我们现在已确认：

• 我们处理的是非平稳时间序列。

• 一阶差分是使时间序列平稳的合适变换技术。

这最终将引导我们进入最后一节。

自相关

一旦我们的数据平稳了，我们可以研究其他关键的时间序列属性：部分自相关和自相关。用正式术语来说：

自相关函数（ACF） 测量时间序列滞后值之间的线性关系。换句话说，它测量时间序列与自身的相关性。[2]

部分自相关函数（PACF） 测量在移除中间相关滞后值的影响后，时间序列中滞后值之间的相关性。这些被称为混杂变量。[3]

这两个指标可以通过称为相关图的统计图进行可视化。但首先，重要的是要对它们有更好的理解。

由于本文主要集中在探索性分析上，并且这些概念是统计预测模型的基础，因此我会简要说明，但请记住，这些都是建立稳固直觉的重要理念，尤其是在处理时间序列时。如果需要全面阅读，我推荐 Kaggle Notebooks 大师Leonie Monigatti的精彩核实文章“时间序列：解读 ACF 和 PACF”。

如上所述，自相关测量时间序列如何与其之前的 q 滞后值相关。你可以将其视为数据子集与其自身向后移动 q 周期的拷贝之间的线性关系的度量。自相关，或 ACF，是确定移动平均（MA）模型阶数 q 的重要指标。

另一方面，部分自相关是时间序列与其 p 滞后版本的相关性，但现在仅涉及其直接影响。例如，如果我想检查 t-3 到 t-1 时间段的部分自相关与当前 t0 值的关系，我不会关注 t-3 如何影响 t-2 和 t-1，或者 t-2 如何影响 t-1。我将专注于 t-3、t-2 和 t-1 对当前时间戳 t0 的直接影响。部分自相关，或 PACF，是确定自回归（AR）模型阶数 p 的重要指标。

理解了这些概念后，我们现在可以回到我们的数据上。由于这两个指标通常一起分析，我们的最后一个函数将会将 PACF 和 ACF 图结合到一个网格图中，这将返回多个变量的相关图。它将利用 statsmodels 的plot_pacf()和plot_acf()函数，并将它们映射到 Matplotlib 的subplots()网格中。

请注意，两个 statsmodels 函数使用了相同的参数，唯一的例外是 plot_pacf() 图的 method 参数。

现在你可以尝试不同的数据聚合，但请记住，当重新采样时间序列时，每个滞后将代表不同的时间回溯。为了说明这一点，我们来分析 2016 年 1 月所有四个站点的 PACF 和 ACF，使用 6 小时聚合的数据集。

图 19. 2016 年 1 月的 PACF 和 ACF 相关图。图像由作者提供。

相关图返回的相关系数范围从 -1.0 到 1.0，并且有一个阴影区域表示显著性阈值。任何超出该范围的值都应被视为统计显著的。

根据上述结果，我们最终可以得出结论，在 6 小时聚合下：

• 滞后 1、2、3（t-6h、t-12h 和 t-18h）以及有时滞后 4（t-24h）具有显著的 PACF。

• 滞后 1 和 4（t-6h 和 t-24h）在大多数情况下显示出显著的 ACF。

并注意一些最终的好实践：

• 避免对长时间段的时间序列进行高粒度的相关图绘制（例如，对具有每小时测量数据集绘制一整年的相关图），因为随着样本量的增加，显著性阈值会缩小到零。

• 我为我们的函数定义了一个 x_label 参数，以便轻松地在 X 轴上标注每个滞后所代表的时间段。通常会看到没有这些信息的相关图，但轻松访问这些信息可以避免对结果的误解。

• Statsmodels 的 plot_acf() 和 plot_pacf() 默认值设置为在图中包含 0 滞后相关系数。由于一个数字与其自身的相关性始终为 1，我已将我们的图从第一个滞后开始绘制，参数为 zero=False。这也改善了 Y 轴的刻度，使我们实际需要分析的滞后更加易读。

这样，我们已经全面探索了我们的时间序列。借助一系列可视化和分析功能，我们可以对数据有一个全面的理解。你还学习了一些在探索时间序列数据集时的最佳实践，以及如何用高质量的图表简洁且精炼地呈现这些数据。

喜欢这个故事吗？

你可以在 Medium 上关注我，获取更多关于数据科学、机器学习、可视化和数据分析的文章。

你还可以在 LinkedIn 和 X 上找到我，我会在这些平台上分享这些内容的简短版本。

## 订阅 Erich Silva 的邮件更新

订阅邮件更新，以便每当 Erich Silva 发布新内容时，你都会收到邮件。如果你还没有 Medium 帐户，将会创建一个…

medium.com

参考

[1] “卫生部门—细颗粒物（PM 2.5）问答。” 访问日期：2022 年 10 月 14 日。 www.health.ny.gov/environmental/indoors/air/pmq_a.htm

[2] Peixeiro, Marco. “3. 随机游走。” 论文。收录于*《Python 中的时间序列预测》*，第 30–58 页。O’Reilly Media，2022 年。

[3] Peixeiro, Marco. “5. 建模自回归过程。” 论文。收录于*《Python 中的时间序列预测》*，第 81–100 页。O’Reilly Media，2022 年。

[4] Peixeiro, Marco. “8. 季节性调整。” 论文。收录于*《Python 中的时间序列预测》*，第 156–179 页。O’Reilly Media，2022 年。

[5] 服务部，公民事务部。“BC 数据目录。” 英属哥伦比亚省。英属哥伦比亚省，2022 年 2 月 2 日。 www2.gov.bc.ca/gov/content/data/bc-data-catalogue

[6] “平稳性和去趋势（ADF/KPSS）。” statsmodels。访问日期：2022 年 10 月 17 日。 www.statsmodels.org/dev/examples/notebooks/generated/stationarity_detrending_adf_kpss.html

计算一组地点坐标的距离矩阵（Python 实现）

原文：towardsdatascience.com/compute-the-distance-matrix-of-a-set-of-sites-from-their-coordinates-in-python-d5fc92a0ba9e?source=collection_archive---------0-----------------------#2023-07-16

轻松估算任意一对地点之间的距离，作为解决一般路由问题的起点

Carlos J. Uribe

·

关注 发表在 Towards Data Science ·13 分钟阅读·2023 年 7 月 16 日

–

由 DALL·E 3 生成的图像，作者的提示：“一张城市网络的地图，每个城市都连接到其他所有城市”

👁️ 这是关于“Python 中的智能决策支持系统”项目的第 4 篇文章， 我鼓励你查看一下，以获取整个项目的一般概述。如果你只对创建距离矩阵感兴趣，这篇文章已经足够，它是自包含的。如果你还想将距离矩阵应用于实际问题，这个系列将对你有兴趣。

本文将继续从sprint 3的最后一个地方开始：为旅行商问题构建优化模型，在给定固定访问地点及其对之间的距离的情况下。在第 4 轮中，我们将暂时从建模中绕开，开发一个具有地理空间功能的类，这将在我们尝试解决一般性旅行商问题时非常有用，即对于任意位置我们可能没有现成的距离数据的问题。我们在上一轮中提出了这个“需求”，将在这一轮中构建一个子系统来满足它。

目录

1. 上一轮迭代回顾

2. 读取输入数据

3. 从位置数据创建距离矩阵

• 3.1. 我是否需要付出额外的努力来获得额外的进展？

• 3.2. 带有的地理定位工具 geopy

• 3.3. 到达要点

• 3.4. 从坐标到距离矩阵

4. 完成！（类内部）

• 4.1. GeoAnalyzer 类设计

• 4.2. 类使用演示

5. 结论（或规划下一轮迭代）

1. 上一轮迭代回顾

在上一篇文章中，即sprint 3，我们进行了概念验证，展示了我们可以解决旅行商问题（TSP）对于一组站点，前提是我们拥有每一对站点之间的距离，作为距离矩阵：

## 用 Python 实现、解决和可视化旅行商问题

学习如何将优化模型从数学翻译到 Python，优化它，并可视化解决方案以快速获得结果…

towardsdatascience.com

我们将距离矩阵视为已给定，因为在那个开发阶段，重点是模型构建，而不是数据获取。但是一旦模型准备好并且在我们的固定地点集上运行良好，我们很快意识到我们需要一种方法来 解决一般的 TSP 问题（任意地点集的问题）。这种泛化是创建真正有用的 MVP 所必需的。因此，我们得出的自然下一步是找到一种自动从我们兴趣点的坐标中获取距离矩阵的方法，这一步我们将在本文中讨论。

这样做，我们的新基本输入将会简单自然得多，只需提供我们想要访问地点的地理坐标：

图 1. 兴趣点的坐标。 （图像由作者提供）

输出将是我们用作 TSP 模型输入的数据框，即输入地点的距离矩阵：

图 2. 给定一组地点的期望距离矩阵。 （图像由作者提供）

为了保持一致性，我们将使用到目前为止考虑的相同巴黎地点。在下一篇文章中，我们将把这个功能与旅行销售员问题的优化模型集成，得到一个更具多功能性的 MVP。

🎯牢记最终目标

让我们稍微回顾一下为什么要做这个。我们期望解决的原始实际问题是我们可以称之为的 旅行游客问题 (TTP)，即为一般游客制定 最佳旅行计划 的问题， 给定她的“个人”数据 (如偏好、预算等) 以及旅行“环境”数据 (如距离、价格、交通方式等)。

由于这个实际问题被认为过于复杂，我们在 第一个冲刺中将其简化为其本质版，以启动解决方案的设计。这个“本质问题”被证明是 旅行销售员问题 (TSP)，在这个问题中，我们将要访问的点视为城市中游客的“兴趣点”。通过本文开发的功能，我们更接近于 TTP 的通用解决方案，以 TSP 作为解决方案的核心。

2. 读取输入数据

我们的基本输入现在是我们旅行中想要访问的地点的地理坐标。我们将“酒店”视为一种不同的地点，因为酒店本身并不是一个“感兴趣的地点”，而是我们必须在多日旅行中停留的地方。我们的酒店选择可能会因不同的旅行或不同的情况而有所不同，而一个城市中的感兴趣的地点则是相对“固定”的地方，许多旅行指南对此意见一致。当我们准备探索更高级的应用时，这种区分的有用性将变得更加明显。

因此，我将我们酒店的坐标存储在一个 CSV 文件 location_hotel.csv 中，将“感兴趣地点”的坐标存储在另一个 CSV 文件 sites_coordinates.csv 中。这两个 CSV 文件具有相同的结构，因此我们将它们读取并合并成一个包含所有地点的数据框：

import pandas as pd

print(f"version pandas: {pd.__version__}")

DATA_FOLDER = ("https://raw.githubusercontent.com/carlosjuribe/"
               "traveling-tourist-problem/main/data")
FILE_LOCATION_HOTEL = "location_hotel.csv"
FILE_LOCATION_SITES = "sites_coordinates.csv"

df_sites = pd.concat([
    pd.read_csv(f"{DATA_FOLDER}/{FILE_LOCATION_SITES}", index_col='site'),
    pd.read_csv(f"{DATA_FOLDER}/{FILE_LOCATION_HOTEL}", index_col='site')
])

display(df_sites)

ℹ️ 如何快速准备自己的位置数据

如果你想使用自己的网站列表来跟随本文，你需要复制我获取坐标的步骤：

1. 前往 Google Maps 并搜索你列表中的每个地点。

2. 每个地点将在地图上显示为一个点。右击每个点。出现的第一个元素是一对数字：你点击的点的纬度和经度。

3. 点击这些数字，它们将被保存到你的剪贴板中，准备粘贴到一个文件中，连同你为该点选择的名称一起。

4. 对所有地点重复步骤 1 到 3，你将得到一个等效于 sites_coordinates.csv 的文件。

这个过程对于小规模的网站集非常有效，但如果你有数百个，甚至几十个地点，它会变得非常繁琐。在 未来的一篇文章 中，我们将创建一种自动化这个手动工作的方式，这叫做地理定位。

3. 从位置数据创建距离矩阵

要构建距离矩阵，我们需要获得任何一对地点之间的距离。这听起来很简单，但“距离”实际上取决于上下文。我们是否考虑由如谷歌地图这样的地图应用程序报告的数字，这些应用程序考虑了街道网络、桥梁、公园，等等？如果是这样，我们考虑步行者行走的距离，还是汽车行驶的距离？或者只是连接这两个点的直线长度？显然，我们有许多可能的距离选择，每种选择的准确度不同。我们必须回答的第一个问题是：在我们的实际问题的特定背景下，我们应该如何定义“距离”，以及在这个阶段？

3.1. 我是否应该多花点时间以获得更多的好处？

很自然地，我们会被诱惑使用准确的数据。最终，我们都知道准确性是内在有价值的，因此我们倾向于追求更准确的数据，越多越好。但我们也必须记住，更准确的数据意味着更复杂的代码和依赖，因此需要更多的开发时间和维护。由于我们遵循敏捷的方法，我们不让最佳成为良好的敌人，因此我们将尽可能简单地开始，然后逐步增加复杂性，仅在必要时。

在需要找到地点之间的距离时，我们可以像很多人一样，直接跳到基于第三方 API 的解决方案，这些方案需要应用密钥、凭据，甚至云服务提供商的信用卡号码。这种方法是可以的，但往往效率不高，因为我们可能会忘记准确的信息带来附加价值，但也会带来附加成本。

👁️ 没有“免费准确性”这种事

记住通常我们总是“付出代价”来获取准确的数据（这与 信息的价值) 紧密相关）是另一个原因，为什么采取敏捷的方法来解决问题是更精简的做法。通过 从简单假设开始 关于“所需的准确度”，并在我们自己的问题数据上验证其有效性 *，我们确保，如果我们最终需要提高数据的准确性，我们将“付出代价”，这种代价是 值得的 (预期的) 改进结果。

所以我们从非常简单的开始。我们有坐标。第一个想法：这些坐标分布在相对于地球半径非常小的地球的地块上，因此我们可以将纬度视为 Y 坐标，将经度视为 X 坐标，在二维平面上进行计算，然后计算欧几里得距离（即通常所说的“直线”）。

• 优点：简单的距离公式，没有新的依赖或数据，地点之间的空间关系被保留。

• 缺点：纬度和经度是无量纲的数字，因此我们在解决问题时得到的数字将不是实际的距离。这意味着我们关心的一些信息，如总旅行距离，即使我们可以获得最佳路线，也不会提供。

缺点胜过优点，因此我们需要一种更复杂的方法（但仍然简单）。第二个想法：将坐标视为它们本身，即地球上的点，但将地球近似为一个球体。一个球体没有我们熟悉的欧几里得几何，因此我们需要一个非平凡的公式，在计算两个点之间的“直线”距离时考虑这种球面几何。所以现在只需使用地球的半径来实现这个公式。我们可以这样做，但我们会依赖一个已经做了这件事的著名库，而且效果更好。

3.2. 使用 **geopy** 的地理位置工具

如果这篇文章系列特别关注地理空间数据科学，那么花时间解释和实现大圆距离的公式将是有价值的，这是一种计算球面上两点“直线”距离的良好基准选项。然而，这篇文章系列关于基于优化的旅游规划系统，因此，我们将依赖于Geopy来完成繁重的工作，而不是自制地理空间工具的公式。这样，我们可以专注于快速找到解决方案。

通过在 Anaconda 提示符（或在我们在第一篇文章中创建的 conda 环境内）运行以下命令来安装它：

conda install -y -c conda-forge geopy=2.3.0

现在，让我们用 geopy 对两个位置进行演示。

3.3. 到达要点

给定两点的坐标，[geodesic](https://geopy.readthedocs.io/en/stable/#geopy.distance.geodesic) 函数可以计算它们之间的地球表面最短距离。在几何学中，测地线 是在给定度量空间上两点之间的最短路径。在我们熟悉的欧几里得空间中，直线就是测地线。在球面空间中，大圆是测地线。Geopy 的 geodesic 函数所考虑的基础“空间”是地球的精确椭球模型。

👁 大圆很棒，但椭圆更棒

*之前我说过，我们将地球视为一个球体，因为这是最简单的可行近似。实际上，地球并不是一个球体，而是一个椭球体，是一种具有更复杂几何形状的固体。现在*geopy*将免于我们为非欧几里得几何编写自己的函数，我们可以升级对地球的近似，并使用更精确的 椭球距离 来代替大圆距离。对于相同的代码行，这确实是免费的精确度，那么为什么不使用呢？

这是一个计算点 1 和点 2 之间椭球距离的函数，单位为米：

from geopy.distance import geodesic

def ellipsoidal_distance(p1, p2) -> float:
    """ Calculate distance (in meters) between p1 and p2, where 
    each point is represented as a tuple (lat, lon) """
    return geodesic(p1, p2).meters

埃菲尔铁塔和卢浮宫之间的距离是多少？

p1 = df_sites.loc['Tour Eiffel']
p2 = df_sites.loc['Louvre']

ellipsoidal_distance(p1, p2)  # output: 3173.119635531859

3173 米，大约 3.2 公里。谷歌地图显示为 3.5 公里。计算的距离比“真实”距离低 8.6%。不过，我们的腿只关心距离的绝对误差，在这种情况下，只比估算距离多走了 330 米，对期望全天在大城市步行的游客来说，这似乎不是一个显著的误差。

那么埃菲尔铁塔和苏弗伦港之间呢？

ellipsoidal_distance(
    df_sites.loc['Tour Eiffel'],
    df_sites.loc['Port de Suffren']
)  # output: 328.3147101635456

328 米，这次比 Google Maps 提供的 350 米低 6%（仅 22 米短）。对于应用一个公式来说，这还不错。正如我们所预期的，点越近，街道出现之类的曲折转弯的机会越小，因此椭球模型产生的误差也越小。对于我们目前的目的来看，显得足够好。

现在我们必须将这个函数应用于所有位置对，从而得到 TSP 模型所需的距离矩阵。

3.4. 从坐标到距离矩阵

这部分很简单，我们只需对所有站点进行两次循环，计算并存储每对之间的距离。下面的函数就是这样做的。注意，距离度量作为可选参数传递，默认是我们之前使用的椭球距离。我们留着以后传递更好的距离度量的可能性。

def compute_distance_matrix(df_sites, dist_metric=ellipsoidal_distance):
    """ Creates an N x N distance matrix from a dataframe of N locations 
    with a latitute column and a longitude column """
    df_dist_matrix = pd.DataFrame(index=df_sites.index, 
                                  columns=df_sites.index)

    for orig, orig_loc in df_sites.iterrows():  # for each origin
        for dest, dest_loc in df_sites.iterrows():  # for each destination
            df_dist_matrix.at[orig, dest] = dist_metric(orig_loc, dest_loc)
    return df_dist_matrix

df_distances = compute_distance_matrix(df_sites)

display(df_distances)

图 3. 使用地球椭球模型得到的距离矩阵。（图像来源：作者）

就这样！正如预期的那样，矩阵的对角线为零，矩阵是对称的。输出数据框的索引和列包含输入站点的名称。

功能演示完毕。现在我们可以更好地方便使用这个函数。让我们以一种便捷的方式将这个功能封装在一个类中，以便于重复使用，更重要的是，为了更容易与我们在上一个冲刺中构建的 TSP 优化模型集成。

4. 总结！(在类内部)

4.1. GeoAnalyzer 类设计

让我们创建一个新的类，GeoAnalyzer，专门用于处理可能出现在路由问题中的地理空间工具。因此，我们的函数compute_distance_matrix自然地作为一个方法嵌入其中。这个类的主要部分目前将包括：

• 包含站点位置的数据框，属性为_df_locations。

• 纯函数ellipsoidal_distance。

• 方法get_distance_matrix，等同于之前的函数compute_distance_matrix，但使用实例属性_df_locations来计算距离。

由于用户可能希望在分析的任何时刻添加新位置，我们包括了一个名为add_locations的方法，该方法接受一个地理坐标的数据框，并将其附加到先前存在的数据框中。

下面可以找到GeoAnalyzer的定义。注意这里还有其他便利的方法和属性未提及。

from typing import Tuple

import pandas as pd
from geopy.distance import geodesic

class GeoAnalyzer:
    """ Utils for geolocation information and processing """  
    _GeoPoint = Tuple[float, float]

    def __init__(self):
        """ Use method `add_locations` to store some locations inside 
        and start using the geo-utilities """
        self._df_locations = pd.DataFrame(columns=['latitude', 'longitude'])

    #####################   distances   #####################
    @staticmethod
    def ellipsoidal_distance(point1: _GeoPoint, point2: _GeoPoint) -> float:
        """ Calculate ellipsoidal distance (in meters) between point1 
        and point2 where each point is represented as a tuple (lat, lon)
        """
        return geodesic(point1, point2).meters
    #########################################################

    @property
    def locations(self):
        return self._df_locations

    @property
    def num_locations(self):
        return len(self._df_locations)

    def add_locations(self, df_locations: pd.DataFrame):
        """ Geo-location data needed for analysis.
        Parameters
        ----------
        df_locations : pd.DataFrame
            Dataframe of geographical coordinates with the first column 
            named 'latitude' and the second column named 'longitude'
        """
        self._name_index = df_locations.index.name
        df_updated = pd.concat([self._df_locations, df_locations.copy()])
        # drop duplicates just in case the user adds repeated locations
        self._df_locations = df_updated.drop_duplicates()

    def get_distance_matrix(self, precision: int = 4) -> pd.DataFrame:
        """ Computes the distance matrix as a dataframe based on the 
        provided location data """
        df_locations = self._df_locations
        dist_metric = self.ellipsoidal_distance  # only distance available

        # initialize matrix df
        df_dist_matrix = pd.DataFrame(index=df_locations.index, 
                                      columns=df_locations.index)
        # for each origin and destination pair, compute distance
        for orig, orig_loc in df_locations.iterrows():
            for dest, dest_loc in df_locations.iterrows():
                distance = round(dist_metric(orig_loc, dest_loc), precision)
                df_dist_matrix.at[orig, dest] = distance

        # a bit of metadata doesn't hurt
        df_dist_matrix.distance_metric = dist_metric.__name__
        df_dist_matrix.index.name = self._name_index
        return df_dist_matrix    

    def __repr__(self):
        """ Display number of currently considered locations """
        return f"{self.__class__.__name__}(n_locs={self.num_locations})"

4.2. 类使用演示

让我们稍微探索一下这个类的主要功能。我们创建一个实例并从巴黎添加我们感兴趣的站点：

geo_analyzer = GeoAnalyzer()
geo_analyzer.add_locations(df_sites)

我们检查此时实例的表示，告知我们已经提供了 9 个位置，我们可以通过属性locations查看详细信息：

display(geo_analyzer)
display(geo_analyzer.locations)

当然，我们可以从对象中提取距离矩阵，到目前为止这已经相当熟悉了：

df_distances = geo_analyzer.get_distance_matrix()

display(df_distances)

最后，如果我们对这些值的来源感到好奇，我们可以从数据框本身进行检查：

print(f"Distance metric used: {df_distances.distance_metric}")
# [Out]: Distance metric used: ellipsoidal_distance

如果有更多的距离度量可用，这将更加有价值，这是我们在未来冲刺中将看到的内容。

5. 结论（或为下一个冲刺做计划）

我们工作的最终结果是一个名为 GeoAnalyzer 的类，具有便捷的方法，帮助我们将旅行推销员问题推广到任意的站点集合。这个推广将是我们下一个冲刺的具体目标，我们将为 TSP 创建一个类似估算器的类，它隐藏了在sprint 2中涉及的模型构建步骤，并以待访问站点的地理坐标作为输入。GeoAnalyzer 类将是这个新估算器类的关键组成部分，实现我们所构建的 TSP 优化模型的真正通用应用。这个新的类似估算器的类，将结合 GeoAnalyzer 和 TSP 模型的通用性，将成为我们解决更一般的旅行游客问题的核心。继续阅读下一个冲刺，了解真正的内容：

## 一种优雅的方法来有效解决旅行推销员问题，使用 Python

以类似 scikit-learn 的方式实现 TSP 模型，以简化路由优化的构建和解决过程…

towardsdatascience.com

随时关注我，向我提问，给我反馈，或通过LinkedIn与我联系。感谢阅读！📈😊

Python 中的并发

原文：towardsdatascience.com/concurrency-in-python-fe8b39edfba5

PYTHON | PROGRAMMING | CONCURRENCY

这是一本关于利用并发执行的力量以及提高 Python 程序性能的初学者指南。

David Farrugia

·发表于 Towards Data Science ·9 分钟阅读·2023 年 5 月 24 日

–

照片由 fabio 提供，来源于 Unsplash

戈登·摩尔在 1965 年做出了一个后来被称为摩尔定律的预测。他指出，微芯片上的晶体管数量每两年将翻一番。此外，摩尔定律还规定，在同一时期内，计算硬件的成本也将减半。

来源：commons.wikimedia.org/wiki/File:Moore%27s_Law_Transistor_Count_1970-2020.png

在今天的技术环境中，计算机设备配备多核 CPU 或多个处理器已很常见。作为开发人员，我们需要编写能够利用这些硬件能力的代码，以提供对用户最优化和高效的解决方案。

什么是并发？

并发是同时执行多个指令序列。

假设我们的系统有一个 2 核心的 CPU。运行非并发代码将导致我们的脚本仅利用一个核心来执行任务，而另一个核心则处于空闲状态。

通过利用并发，我们可以在两个核心上同时执行任务，从而提高性能并减少等待时间。

然而，并发的一个缺点是我们无法保证任务执行的顺序。

因此，重要的是各种任务的执行顺序不应影响结果，并且任务应尽可能共享较少的资源。

为了缓解这一点，需要协调共享资源，这会增加复杂性。

共享资源越多，复杂性越高。

并发的类型

并行编程

并行编程是将主要任务分解为较小的子任务的实践。

这些子任务然后被分配到不同的线程或进程中，并在多个核心上同时执行。

相反，在单核编程中，只使用一个核心来执行任务，其他核心则处于空闲状态或可供其他程序或任务使用。

这可能导致资源利用效率降低。

通过利用并行编程，可以使用所有核心来提高整体性能。

并行编程的示例。来源：作者

由于并行编程利用了多个核心，因此它特别适用于 CPU 绑定的任务，也称为 CPU 密集型任务。这些任务只能通过增加更多的处理器来加速。此类任务的例子包括搜索算法和数学运算。

异步编程

在异步编程中，也称为多线程，主线程将子任务发送给一个演员，这可以是另一个线程或设备。

主线程随后继续执行其他工作，而不是等待响应。当子任务完成时，演员通知主线程并触发回调函数来处理结果。

在 Python 中，我们使用一个称为“future”的对象来代替回调函数，该对象表示尚未完成操作的结果。

异步编程的示例。来源：作者

根据程序的结构，主线程将等待子任务完成或在稍后时间检查。

异步编程特别适用于 I/O 绑定任务。I/O 绑定任务是 CPU 无法单独执行的操作，而是依赖于磁盘访问或网络等外部因素。I/O 绑定任务的例子包括从/向存储读/写和发出 API 请求。

并行与异步编程

在任何给定的时间点：

并行编程 — 更快地完成任务

异步编程 — 做更多的事情

Python 中的多线程

一个常用的例子来说明多线程的好处是同时下载多个网络上的图像。我们可以利用这个例子来理解多线程的有效性以及如何在 Python 中实现它。

import os
import time
from urllib.parse import urlparse
from urllib.request import urlretrieve
from typing import List
from numpy import round

IMGS_URL_LIST = \
    ['https://dl.dropboxusercontent.com/s/2fu69d8lfesbhru/pexels-photo-48603.jpeg',
     'https://dl.dropboxusercontent.com/s/zch88m6sb8a7bm1/pexels-photo-134392.jpeg',
     'https://dl.dropboxusercontent.com/s/lsr6dxw5m2ep5qt/pexels-photo-135130.jpeg',
     'https://dl.dropboxusercontent.com/s/6xinfm0lcnbirb9/pexels-photo-167300.jpeg',
     'https://dl.dropboxusercontent.com/s/2dp2hli32h9p0y6/pexels-photo-167921.jpeg',
     'https://dl.dropboxusercontent.com/s/fjb1m3grcrceqo2/pexels-photo-173125.jpeg',
     'https://dl.dropboxusercontent.com/s/56u8p4oplagc4bp/pexels-photo-185934.jpeg',
     'https://dl.dropboxusercontent.com/s/2s1x7wz4sdvxssr/pexels-photo-192454.jpeg',
     'https://dl.dropboxusercontent.com/s/1gjphqnllzm10hh/pexels-photo-193038.jpeg',
     'https://dl.dropboxusercontent.com/s/pcjz40c8pxpy057/pexels-photo-193043.jpeg',
     'https://dl.dropboxusercontent.com/s/hokdfk7y8zmwe96/pexels-photo-207962.jpeg',
     'https://dl.dropboxusercontent.com/s/k2tk2co7r18juy7/pexels-photo-247917.jpeg',
     'https://dl.dropboxusercontent.com/s/m4xjekvqk4rksbx/pexels-photo-247932.jpeg',
     'https://dl.dropboxusercontent.com/s/znmswtwhcdbpc10/pexels-photo-265186.jpeg',
     'https://dl.dropboxusercontent.com/s/jgb6n4esquhh4gu/pexels-photo-302899.jpeg',
     'https://dl.dropboxusercontent.com/s/rjuggi2ubc1b3bk/pexels-photo-317156.jpeg',
     'https://dl.dropboxusercontent.com/s/cpaog2nwplilrz9/pexels-photo-317383.jpeg',
     'https://dl.dropboxusercontent.com/s/16x2b6ruk18gji5/pexels-photo-320007.jpeg',
     'https://dl.dropboxusercontent.com/s/xqzqzjkcwl52en0/pexels-photo-322207.jpeg',
     'https://dl.dropboxusercontent.com/s/frclthpd7t8exma/pexels-photo-323503.jpeg',
     'https://dl.dropboxusercontent.com/s/7ixez07vnc3jeyg/pexels-photo-324030.jpeg',
     'https://dl.dropboxusercontent.com/s/1xlgrfy861nyhox/pexels-photo-324655.jpeg',
     'https://dl.dropboxusercontent.com/s/v1b03d940lop05d/pexels-photo-324658.jpeg',
     'https://dl.dropboxusercontent.com/s/ehrm5clkucbhvi4/pexels-photo-325520.jpeg',
     'https://dl.dropboxusercontent.com/s/l7ga4ea98hfl49b/pexels-photo-333529.jpeg',
     'https://dl.dropboxusercontent.com/s/rleff9tx000k19j/pexels-photo-341520.jpeg'
    ]

def download_images(img_url_list: List[str]) -> None:
    # validate inputs
    if not img_url_list:
        return
    os.makedirs('images', exist_ok=True)

    # get time in seconds
    start = time.perf_counter()

    # for every url in our list, we parse the url and download its contents
    for img_num, url in enumerate(img_url_list):
        urlretrieve(url, f'images{os.path.sep}{img_num+1}')

    print(f"Retrieved {len(img_url_list)} images took {round(time.perf_counter() - start, 2)} seconds")

download_images(IMGS_URL_LIST)

在上述单线程脚本中，我们设置了一个函数（download_images），用于检索公开托管在 Dropbox 上的一些图像。

这个脚本下载 26 张图像花费了 22.06 秒。

为了改善这一点，我们可以转变编程思路，使用不同的方法来构建一个利用并发的脚本。

我们可以将逻辑分成两个主要函数：目标函数和运行函数，而不是编写一个循环遍历每个 URL 并检索其内容的函数。

目标函数的作用是封装处理单个 URL 所需的逻辑。

由于我们希望为每个 URL 拥有一个线程，因此我们需要为每个线程提供处理 URL 的知识。

然后，运行函数用于触发每个 URL 的新线程并存储它们的结果。

在运行函数中，我们需要指示主线程为每个 URL 创建并启动一个子线程。

我们可以通过遍历传递的 URL 并在 Python 中创建和启动一个新线程来做到这一点，如下所示：

from threading import Thread

t = Thread(target=<TARGET_FUNCTION>, args=(<SOME_ARGS>))
t.start()

在我们的例子中，我们希望等待所有图像下载完成后再继续程序。为此，我们可以指示主线程等待所有子线程完成，然后再继续程序执行，通过调用 join()函数。此函数将子线程与主线程连接在一起，主线程在所有连接的线程执行完成之前不会继续。

from threading import Thread

# store our threads
threads = []

t = Thread(target=<TARGET_FUNCTION>, args=(<SOME_ARGS>))
t.start()
threads.append(t)

# join all child threads to the main thread
for thread in threads:
  thread.join()

通过使用 join()函数，我们可以确保程序在继续执行脚本的下一步之前等待所有子线程完成。

将原始脚本调整为使用多线程大致如下：

import threading
import os
import time
from urllib.parse import urlparse
from urllib.request import urlretrieve
from typing import List
from numpy import round

IMGS_URL_LIST = \
    ['https://dl.dropboxusercontent.com/s/2fu69d8lfesbhru/pexels-photo-48603.jpeg',
     'https://dl.dropboxusercontent.com/s/zch88m6sb8a7bm1/pexels-photo-134392.jpeg',
     'https://dl.dropboxusercontent.com/s/lsr6dxw5m2ep5qt/pexels-photo-135130.jpeg',
     'https://dl.dropboxusercontent.com/s/6xinfm0lcnbirb9/pexels-photo-167300.jpeg',
     'https://dl.dropboxusercontent.com/s/2dp2hli32h9p0y6/pexels-photo-167921.jpeg',
     'https://dl.dropboxusercontent.com/s/fjb1m3grcrceqo2/pexels-photo-173125.jpeg',
     'https://dl.dropboxusercontent.com/s/56u8p4oplagc4bp/pexels-photo-185934.jpeg',
     'https://dl.dropboxusercontent.com/s/2s1x7wz4sdvxssr/pexels-photo-192454.jpeg',
     'https://dl.dropboxusercontent.com/s/1gjphqnllzm10hh/pexels-photo-193038.jpeg',
     'https://dl.dropboxusercontent.com/s/pcjz40c8pxpy057/pexels-photo-193043.jpeg',
     'https://dl.dropboxusercontent.com/s/hokdfk7y8zmwe96/pexels-photo-207962.jpeg',
     'https://dl.dropboxusercontent.com/s/k2tk2co7r18juy7/pexels-photo-247917.jpeg',
     'https://dl.dropboxusercontent.com/s/m4xjekvqk4rksbx/pexels-photo-247932.jpeg',
     'https://dl.dropboxusercontent.com/s/znmswtwhcdbpc10/pexels-photo-265186.jpeg',
     'https://dl.dropboxusercontent.com/s/jgb6n4esquhh4gu/pexels-photo-302899.jpeg',
     'https://dl.dropboxusercontent.com/s/rjuggi2ubc1b3bk/pexels-photo-317156.jpeg',
     'https://dl.dropboxusercontent.com/s/cpaog2nwplilrz9/pexels-photo-317383.jpeg',
     'https://dl.dropboxusercontent.com/s/16x2b6ruk18gji5/pexels-photo-320007.jpeg',
     'https://dl.dropboxusercontent.com/s/xqzqzjkcwl52en0/pexels-photo-322207.jpeg',
     'https://dl.dropboxusercontent.com/s/frclthpd7t8exma/pexels-photo-323503.jpeg',
     'https://dl.dropboxusercontent.com/s/7ixez07vnc3jeyg/pexels-photo-324030.jpeg',
     'https://dl.dropboxusercontent.com/s/1xlgrfy861nyhox/pexels-photo-324655.jpeg',
     'https://dl.dropboxusercontent.com/s/v1b03d940lop05d/pexels-photo-324658.jpeg',
     'https://dl.dropboxusercontent.com/s/ehrm5clkucbhvi4/pexels-photo-325520.jpeg',
     'https://dl.dropboxusercontent.com/s/l7ga4ea98hfl49b/pexels-photo-333529.jpeg',
     'https://dl.dropboxusercontent.com/s/rleff9tx000k19j/pexels-photo-341520.jpeg'
    ]

# this is our target function
def download_image(url: str, img_num: int) -> None:
    urlretrieve(url, f'images{os.path.sep}{img_num+1}')

def download_images(img_url_list: List[str]) -> None:
    # validate inputs
    if not img_url_list:
        return
    os.makedirs('images', exist_ok=True)

    # get time in seconds
    start = time.perf_counter()

    # create a list to store all of our threads
    threads = []

    # for every url in our list, we parse the url and download its contents
    for img_num, url in enumerate(img_url_list):
        # create a new thread
        t = threading.Thread(target=download_image, args=(url, img_num))

        # start the new thread
        t.start()

        # add the new thread to our list of threads
        threads.append(t)

    # here we instruct the main thread to wait for all child threads to complete before proceeding
    for thread in threads:
        thread.join()

print(f"Retrieved {len(img_url_list)} images took {round(time.perf_counter() - start, 2)} seconds")

download_images(IMGS_URL_LIST)

这个脚本仅需 2.92 秒即可检索相同的 26 张图像！这只是使用单线程代码所需时间的 13.24%。

令人印象深刻，不是吗？

理解执行流程

在深入了解执行流程之前，我们必须首先了解不同的线程状态。

不同的线程状态。来源：作者

• 新建 — 一个新创建的线程

• 就绪 — 线程准备好执行

• 运行中 — 线程正在执行

• 阻塞 — 线程被阻塞执行（等待资源变得可用）

• 完成 — 线程执行完毕

当我们开始执行时，只有一个线程可用：主线程。

主线程开始执行代码，直到它遇到生成新线程的指令。这时，我们现在有两个线程：主线程和一个子线程。

子线程从“新建”状态移动到“就绪”状态。当主线程执行 start()函数时，子线程移动到“运行中”状态并开始执行。

主线程变为“阻塞”状态，直到子线程处于“完成”状态（由于join()函数）。

在更复杂的任务中，尤其是当线程访问共享资源时，控制哪个线程在任何给定时间可以访问哪个资源是很重要的。

这个过程通过使用信号量、锁和其他变量来实现。但我们可以在其他时间深入探讨这个话题。

Python 中的多处理

在 Python 中创建新进程很简单，类似于使用线程。

我们使用 multiprocessing 包中的 Process 函数创建一个新进程。

我们为 Process 对象指定目标函数并传递相应的参数。

start 函数开始运行进程，我们可以使用 join 函数来指示主进程在继续到下一个代码块之前等待所有子进程完成。

与线程不同，进程可以通过调用该进程的 terminate()函数来终止。

如果进程正在运行，is_alive()函数将返回 True，如果不在运行则返回 False。

然而，终止进程可能会导致任何共享资源处于不一致状态，使其对其他进程不可用。因此，在终止进程时要小心。

并行化程序的主要优势是能够利用所有可用的 CPU 来执行我们的任务。

进程池允许我们以几种方式将任务分配给一组称为工作者的进程。

初始化进程池时，我们需要指定池的大小（即工作线程的数量），并且可以选择指定一个初始化函数。

初始化函数将在每个新创建的进程中执行。

选择工作线程数量的一个好起点是默认使用系统可用的核心数量。

我们使用 map 函数在一个可迭代对象上执行目标函数。

 import multiprocessing

# create an initialiser function which prints the name of the process created
def start_process() -> None:
    print(f'Process {multiprocessing.current_process().name} started!')

# create a sample target function which squares a passed integer
def target_function(num: int) -> int:
    return num**2

# get the number of CPUs available    
pool_size = multiprocessing.cpu_count()

# initialise our process pool
pool = multiprocessing.Pool(processes=pool_size, initializer=start_process)

# map target function to our iterable
result = pool.map(target_function, range(200))

# when finished, close the pool
pool.close()

# wait for child processes to complete
pool.join()

# access our results
print(result)

在这种情况下，任务如此简单，以至于多进程会导致比单处理器流程更慢的执行速度。

慢速主要是由于与进程间通信（IPC）相关的开销，包括创建和启动各个进程所花费的时间。

因此，使用 Pool 功能的一些指示包括高 CPU 要求和需要遍历的长可迭代对象。

结论

并发是一种编程范式，它允许我们的程序同时做更多的事情。

在这篇文章中，我们讨论了不同类型的并发，并提供了示例来帮助你开始使用 Python 进行并发。

建议在日常脚本中尝试这些示例，以实践并欣赏你的程序变得更加高效的方式。

当涉及到并发时，还有许多其他高级概念可以探索，但那是另一个话题！

你喜欢这篇文章吗？只需$5/月，你就可以成为会员，解锁对 Medium 的无限访问。你将直接支持我和你所有其他喜欢的 Medium 作者。非常感谢！

## 加入 Medium，使用我的推荐链接 - David Farrugia

获取对我所有⚡高级⚡内容的独占访问权限，并在 Medium 上无限制地访问所有内容。通过购买我…

david-farrugia.medium.com/membership?source=post_page-----fe8b39edfba5--------------------------------

想要联系我？

我很想听听你对这个话题的看法，或者任何关于 AI 和数据的想法。

如果你希望联系我，可以发邮件至davidfarrugia53@gmail.com。

LinkedIn