scanpy源码浅析：或谈python面向对象结构

作者: 周运来就是我 | 来源:发表于2020-04-02 00:03 被阅读0次

我们已经知道scanpy是一个由Python组织的功能强大的单细胞数据分析工具，用的顺手的同时，我们不禁要问：她是如何组织的？

于是，我们尝试阅读他的源代码。在此之前我们要知道scanpy的源码在哪里放着：

python 
Python 3.7.4 (default, Aug 13 2019, 20:35:49) 
[GCC 7.3.0] :: Anaconda, Inc. on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> import scanpy as sc 
>>> sc
<module 'scanpy' from 'pathto/lib/python3.7/site-packages/scanpy/__init__.py'>

映入我们眼帘的是这样的结构，第一步就是要理解这里的文件都是什么作用及其调用关系：

.
├── api【调用接口】
├── cli.py
├── _compat.py
├── datasets 【数据集】
├── external【外部接口】
├── get.py
├── __init__.py
├── logging.py
├── __main__.py
├── neighbors【计算图结构】
├── plotting【绘图模块】
├── preprocessing【预处理模块】
├── __pycache__
├── queries【查询，富集等】
├── readwrite.py
├── _settings.py
├── sim_models
├── tools
└── _utils.py

首先我们要理解一下init.py这个文件。在scanpy中我们看到主要的面向用户的函数都在这里的，这里的 __init__.py是这样的：

# some technical stuff
import sys
from ._utils import pkg_version, check_versions, annotate_doc_types

__author__ = ', '.join([
    'Alex Wolf',
    'Philipp Angerer',
    'Fidel Ramirez',
    'Isaac Virshup',
    'Sergei Rybakov',
    'Gokcen Eraslan',
    'Tom White',
    'Malte Luecken',
    'Davide Cittaro',
    'Tobias Callies',
    'Marius Lange',
    'Andrés R. Muñoz-Rojas',
])
__email__ = ', '.join([
    'f.alex.wolf@gmx.de',
    'philipp.angerer@helmholtz-muenchen.de',
    # We don’t need all, the main authors are sufficient.
])
try:
    from setuptools_scm import get_version
    __version__ = get_version(root='..', relative_to=__file__)
    del get_version
except (LookupError, ImportError):
    __version__ = str(pkg_version(__name__))

check_versions()
del pkg_version, check_versions

# the actual API
from ._settings import settings, Verbosity  # start with settings as several tools are using it
from . import tools as tl
from . import preprocessing as pp
from . import plotting as pl
from . import datasets, logging, queries, external, get

from anndata import AnnData
from anndata import read_h5ad, read_csv, read_excel, read_hdf, read_loom, read_mtx, read_text, read_umi_tools
from .readwrite import read, read_10x_h5, read_10x_mtx, write, read_visium
from .neighbors import Neighbors

set_figure_params = settings.set_figure_params

# has to be done at the end, after everything has been imported
annotate_doc_types(sys.modules[__name__], 'scanpy')
del sys, annotate_doc_types

__init__.py该文件的作用就是相当于把自身整个文件夹当作一个包来管理，每当有外部import的时候，就会自动执行里面的函数。通常init.py 文件为空，但是我们还可以为它增加其他的功能。我们在导入一个包时，实际上是导入了它的init.py文件。这样我们可以在init.py文件中批量导入我们所需要的模块，不再需要一一导入。

api

整个scanpy中，api是其他模块的接口，也就是api中包含了其他模块的调用方式。

tree api/
api/
├── datasets.py
├── export_to.py
├── __init__.py
├── logging.py
├── pl.py
├── pp.py
├── __pycache__
│   ├── datasets.cpython-37.pyc
│   ├── export_to.cpython-37.pyc
│   ├── __init__.cpython-37.pyc
│   ├── logging.cpython-37.pyc
│   ├── pl.cpython-37.pyc
│   ├── pp.cpython-37.pyc
│   ├── queries.cpython-37.pyc
│   └── tl.cpython-37.pyc
├── queries.py
└── tl.py

1 directory, 16 files

api就像一个管理器，没有明显地定义具体的功能。具体的方法在下面的模块中定义，主要有三个核心功能：

pp 【预处理】
tl 【计算】
pl【绘图】

以及datasets、logging、queries、export_to四个不常用的功能。

datasets/
├── _ebi_expression_atlas.py
├── __init__.py
└── __pycache__
    ├── _ebi_expression_atlas.cpython-37.pyc
    └── __init__.cpython-37.pyc

1 directory, 4 files

datasets中__init__.py是这样的，调用不同的数据，这些数据不在本地第一次使用函数时会下载。

"""Builtin Datasets.
"""
from ._datasets import (
    blobs,
    burczynski06,
    krumsiek11,
    moignard15,
    paul15,
    toggleswitch,
    pbmc68k_reduced,
    pbmc3k,
    pbmc3k_processed,
    visium_sge,
)
from ._ebi_expression_atlas import ebi_expression_atlas

preprocessing中定义了预处理的相关模块：highly_variable_genes、normalization、qc等。

tree   preprocessing/ 
preprocessing/
├── _combat.py
├── _deprecated
│   ├── highly_variable_genes.py
│   ├── __init__.py
│   └── __pycache__
│       ├── highly_variable_genes.cpython-37.pyc
│       └── __init__.cpython-37.pyc
├── _distributed.py
├── _docs.py
├── _highly_variable_genes.py
├── __init__.py
├── _normalization.py
├── __pycache__
│   ├── _combat.cpython-37.pyc
│   ├── _distributed.cpython-37.pyc
│   ├── _docs.cpython-37.pyc
│   ├── _highly_variable_genes.cpython-37.pyc
│   ├── __init__.cpython-37.pyc
│   ├── _normalization.cpython-37.pyc
│   ├── _qc.cpython-37.pyc
│   ├── _recipes.cpython-37.pyc
│   ├── _simple.cpython-37.pyc
│   ├── _utils.cpython-37.pyc
│   ├── _utils.sparse_mean_var_minor_axis-49.py37m.1.nbc
│   └── _utils.sparse_mean_var_minor_axis-49.py37m.nbi
├── _qc.py
├── _recipes.py
├── _simple.py
└── _utils.py

如均一化脚本：

from typing import Optional, Union, Iterable, Dict

import numpy as np
from anndata import AnnData
from scipy.sparse import issparse
from sklearn.utils import sparsefuncs

from .. import logging as logg
from .._compat import Literal


def _normalize_data(X, counts, after=None, copy=False):
    X = X.copy() if copy else X
    if issubclass(X.dtype.type, (int, np.integer)):
        X = X.astype(np.float32)  # TODO: Check if float64 should be used
    counts = np.asarray(counts)  # dask doesn't do medians
    after = np.median(counts[counts>0], axis=0) if after is None else after
    counts += (counts == 0)
    counts = counts / after
    if issparse(X):
        sparsefuncs.inplace_row_scale(X, 1/counts)
    else:
        np.divide(X, counts[:, None], out=X)
    return X

在scanpy函数中，我们经常看到，这样的定义结构：

def read(
    filename: Union[Path, str],
    backed: Optional[Literal['r', 'r+']] = None,
    sheet: Optional[str] = None,
    ext: Optional[str] = None,
    delimiter: Optional[str] = None,
    first_column_names: bool = False,
    backup_url: Optional[str] = None,
    cache: bool = False,
    cache_compression: Union[Literal['gzip', 'lzf'], None, Empty] = _empty,
    **kwargs,
) -> AnnData:
    """\
    Read file and return :class:`~anndata.AnnData` object.

    To speed up reading, consider passing ``cache=True``, which creates an hdf5
    cache file.

    Parameters
    ----------
    filename
        If the filename has no file extension, it is interpreted as a key for
        generating a filename via ``sc.settings.writedir / (filename +
        sc.settings.file_format_data)``.  This is the same behavior as in
        ``sc.read(filename, ...)``.
    backed

def A() ->B:"""C"""。这是什么意思呢？

在一些Python的工程项目中，我们会看到函数参数中会有冒号，有的函数后面会跟着一个箭头，你可能会疑惑，这些都是什么东西？

其实函数参数中的冒号是参数的类型建议符，告诉程序员希望传入的实参的类型。函数后面跟着的箭头是函数返回值的类型建议符，用来说明该函数返回的值是什么类型。

更官方的解释：此为type hints，是Python 3.5新加的功能，作用如上所述，官方文档为 https://www.python.org/dev/peps/pep-0484/

值得注意的是，类型建议符并非强制规定和检查，也就是说即使传入的实际参数与建议参数不符，也不会报错。我认为类型建议符的作用更多的体现在软件工程方面：在多人合作的时候，我们对他人开发的代码并不熟悉，没有对类型的解释说明的话，往往需要花费更多的时间才能看出函数的参数和返回值是什么类型，有了说明符，可以方便程序员理解函数的输入与输出（具体涉及到的工作，比如静态分析与代码重构）。

Python函数参数中的冒号与箭头

在tools中写了scanpy的核心计算函数，比如paga、louvain、score_genes、tsne、umap等。

tree  tools/  
tools/
├── _dendrogram.py
├── _diffmap.py
├── _dpt.py
├── _draw_graph.py
├── _embedding_density.py
├── _ingest.py
├── __init__.py
├── _leiden.py
├── _louvain.py
├── _marker_gene_overlap.py
├── _paga.py
├── _pca.py
├── __pycache__
│   ├── _dendrogram.cpython-37.pyc
│   ├── _diffmap.cpython-37.pyc
│   ├── _dpt.cpython-37.pyc
│   ├── _draw_graph.cpython-37.pyc
│   ├── _embedding_density.cpython-37.pyc
│   ├── _ingest.cpython-37.pyc
│   ├── __init__.cpython-37.pyc
│   ├── _leiden.cpython-37.pyc
│   ├── _louvain.cpython-37.pyc
│   ├── _marker_gene_overlap.cpython-37.pyc
│   ├── _paga.cpython-37.pyc
│   ├── _pca.cpython-37.pyc
│   ├── _rank_genes_groups.cpython-37.pyc
│   ├── _score_genes.cpython-37.pyc
│   ├── _sim.cpython-37.pyc
│   ├── _top_genes.cpython-37.pyc
│   ├── _tsne.cpython-37.pyc
│   ├── _tsne_fix.cpython-37.pyc
│   ├── _umap.cpython-37.pyc
│   ├── _utils_clustering.cpython-37.pyc
│   └── _utils.cpython-37.pyc
├── _rank_genes_groups.py
├── _score_genes.py
├── _sim.py
├── _top_genes.py
├── _tsne_fix.py
├── _tsne.py
├── _umap.py
├── _utils_clustering.py
└── _utils.py

1 directory, 42 files

如umap的计算：

# 这里省略210行代码。

。。。代码去谈恋爱了。

趁着代码谈恋爱的阶段我们来介绍一下python中的_(单下划线)以及__（双下划线）。

单下划线:单下划线的变量是一种程序员之间美丽的约定——只要是这种变量就不要随便在类外部去访问它！！！
但是如果我们在导入模块时来看这个单下划线开头的变量，那就不一样了，在这里这种特殊名字的变量就变成了类似一种某个模块的“私有”变量，因为我们在使用from 模块名 import *语句导入模块时，这些单下划线开头的变量默认是不会被导入的，所以实际上这个单下划线对python的解释器有了影响。
双下划线：解析器自动转换为：_类名_成员名，代替原有成员，访问需要在原有成员名字前加上类名。如:Python自动将__name 解释成\ _student__name,我们可以用\ _student__name访问.
python中的单下划线,双下划线以及两端双下划线

`_ 的含义`

在python的类中没有真正的私有属性或方法，没有真正的私有化。

但为了编程的需要，我们常常需要区分私有方法和共有方法以方便管理和调用。那么在Python中如何做呢？

一般Python约定加了下划线 _ 的属性和方法为私有方法或属性，以提示该属性和方法不应在外部调用，也不会被from ModuleA import * 导入。如果真的调用了也不会出错,但不符合规范。

`__ 的含义`

Python中的__和一项称为name mangling的技术有关，name mangling (又叫name decoration命名修饰).在很多现代编程语言中,这一技术用来解决需要唯一名称而引起的问题,比如命名冲突/重载等.

Python中双下划线开头,是为了不让子类重写该属性方法.通过类的实例化时自动转换,在类中的双下划线开头的属性方法前加上”_类名”实现.

构建图结构

tree neighbors/
neighbors/
├── __init__.py
└── __pycache__
    └── __init__.cpython-37.pyc

1 directory, 2 files

绘图函数

tree  plotting/  
plotting/
├── _anndata.py
├── _docs.py
├── __init__.py
├── palettes.py
├── _preprocessing.py
├── __pycache__
│   ├── _anndata.cpython-37.pyc
│   ├── _docs.cpython-37.pyc
│   ├── __init__.cpython-37.pyc
│   ├── palettes.cpython-37.pyc
│   ├── _preprocessing.cpython-37.pyc
│   ├── _qc.cpython-37.pyc
│   ├── _rcmod.cpython-37.pyc
│   └── _utils.cpython-37.pyc
├── _qc.py
├── _rcmod.py
├── _tools
│   ├── __init__.py
│   ├── paga.py
│   ├── __pycache__
│   │   ├── __init__.cpython-37.pyc
│   │   ├── paga.cpython-37.pyc
│   │   └── scatterplots.cpython-37.pyc
│   └── scatterplots.py
└── _utils.py

3 directories, 22 files

我们看看著名的stacked_violin函数：

@_doc_params(show_save_ax=doc_show_save_ax, common_plot_args=doc_common_plot_args)
def stacked_violin(
    adata: AnnData,
    var_names: Union[_VarNames, Mapping[str, _VarNames]],
    groupby: Optional[str] = None,
    log: bool = False,
    use_raw: Optional[bool] = None,
    num_categories: int = 7,
    figsize: Optional[Tuple[float, float]] = None,
    dendrogram: Union[bool, str] = False,
    gene_symbols: Optional[str] = None,
    var_group_positions: Optional[Sequence[Tuple[int, int]]] = None,
    var_group_labels: Optional[Sequence[str]] = None,
    standard_scale: Optional[Literal['var', 'obs']] = None,
    var_group_rotation: Optional[float] = None,
    layer: Optional[str] = None,
    stripplot: bool = False,
    jitter: Union[float, bool] = False,
    size: int = 1,
    scale: Literal['area', 'count', 'width'] = 'width',
    order: Optional[Sequence[str]] = None,
    swap_axes: bool = False,
    show: Optional[bool] = None,
    save: Union[bool, str, None] = None,
    row_palette: str = 'muted',
    ax: Optional[_AxesSubplot] = None,
    **kwds,
):

这里我们插播一条python的知识点：装饰器。

装饰器就是一个函数对另一个函数的管理和持有。

装饰器背后的主要动机源自python面向对象编程，装饰器是在函数调用之上的修饰，这些修饰仅是当声明一个函数或者方法的时候，才会应用的额外调用。装饰器的语法以@开头，接着是装饰器显式的名字和可选的参数。紧跟着装饰器声明的是被修饰的函数，和修饰函数的可选参数。

装饰器实际上是一个函数，他们接受函数对象，但他们是怎么处理那些函数？

当包装一个函数的时候，你最终会调用它，最棒的是我们能够在包装的环境下在合适的时机调用它，我们在执行函数之前，可以运行那些预备代码，也可以在执行代码之后做个清理工作。所以一个装饰器函数，很可能在里面是这样一些代码：它定义了某个函数并在定义内的某处嵌入了对目标函数的调用或者至少一些引用。从本质上看，这些特征引入了java开发者称之为AOP的概念，可以考虑在装饰器中置入通用功能的代码来降低程序复杂度。

可以用装饰器来：

引入日志
增加计时逻辑来检测性能
给函数加入事物能力。

查询功能模块

tree   queries/  
queries/
├── __init__.py
├── __pycache__
│   ├── __init__.cpython-37.pyc
│   └── _queries.cpython-37.pyc
└── _queries.py

1 directory, 4 files

from ._queries import (
    biomart_annotations,
    gene_coordinates,
    mitochondrial_genes,
)  # Biomart queries
from ._queries import enrich  # gprofiler queries

scanpy也是可以做富集分析的啊。我去，python能做的都可以在python中做。

sc.queries.enrich(['Klf4', 'Pax5', 'Sox2', 'Nanog'], org="hsapiens")
Out[166]: 
   source  ...                                            parents
0    REAC  ...                               [REAC:R-HSA-1266738]
1   GO:BP  ...                           [GO:0001704, GO:0045165]
2    REAC  ...                                [REAC:R-HSA-452723]
3   GO:MF  ...                           [GO:0000981, GO:0001216]
4   GO:MF  ...                                       [GO:0003700]
5      WP  ...                                        [WP:000000]
6   GO:BP  ...                                       [GO:0019827]
7      WP  ...                                        [WP:000000]
8   GO:BP  ...                           [GO:0001708, GO:0001711]
9   GO:BP  ...                           [GO:0009653, GO:0009790]

external 外部功能模块

在这里，我们可以发挥想象力了。这里可以成为我们自己的接口，接入我们自己的分析模块。

tree  external/ 
external/
├── exporting.py
├── __init__.py
├── pl.py
├── pp
│   ├── _bbknn.py
│   ├── _dca.py
│   ├── __init__.py
│   ├── _magic.py
│   ├── _mnn_correct.py
│   └── __pycache__
│       ├── _bbknn.cpython-37.pyc
│       ├── _dca.cpython-37.pyc
│       ├── __init__.cpython-37.pyc
│       ├── _magic.cpython-37.pyc
│       └── _mnn_correct.cpython-37.pyc
├── __pycache__
│   ├── exporting.cpython-37.pyc
│   ├── __init__.cpython-37.pyc
│   └── pl.cpython-37.pyc
└── tl
    ├── __init__.py
    ├── _palantir.py
    ├── _phate.py
    ├── _phenograph.py
    ├── __pycache__
    │   ├── __init__.cpython-37.pyc
    │   ├── _palantir.cpython-37.pyc
    │   ├── _phate.cpython-37.pyc
    │   ├── _phenograph.cpython-37.pyc
    │   ├── _pypairs.cpython-37.pyc
    │   └── _trimap.cpython-37.pyc
    ├── _pypairs.py
    └── _trimap.py

5 directories, 28 files

哇，许多可以分析单细胞的工具啊。

from ._pypairs import cyclone, sandbag
from ._phate import phate
from ._phenograph import phenograph
from ._palantir import palantir
from ._trimap import trimap
from ._harmony_timeseries import harmony_timeseries
from ._sam import sam
from ._wishbone import wishbone

sim_models

tree  sim_models/  
sim_models/
├── __init__.py
└── __pycache__
    └── __init__.cpython-37.pyc

1 directory, 2 files

pycache

tree __pycache__/
__pycache__/
├── cli.cpython-37.pyc
├── _compat.cpython-37.pyc
├── get.cpython-37.pyc
├── __init__.cpython-37.pyc
├── logging.cpython-37.pyc
├── __main__.cpython-37.pyc
├── readwrite.cpython-37.pyc
├── _settings.cpython-37.pyc
└── _utils.cpython-37.pyc

0 directories, 9 files

这里再插播一条python的知识点：__pycache__的作用。

Python程序运行时不需要编译成二进制代码，而直接从源码运行程序，简单来说是，Python解释器将源码转换为字节码，然后再由解释器来执行这些字节码。

解释器的具体工作：

1、完成模块的加载和链接；
2、将源代码编译为PyCodeObject对象(即字节码)，写入内存中，供CPU读取；
3、从内存中读取并执行，结束后将PyCodeObject写回硬盘当中，也就是复制到.pyc或.pyo文件中，以保存当前目录下所有脚本的字节码文件。

之后若再次执行该脚本，它先检查【本地是否有上述字节码文件】和【该字节码文件的修改时间是否在其源文件之后】，是就直接执行，否则重复上述步骤。

那么，__pycache__文件夹的意义何在呢？
因为第一次执行代码的时候，Python解释器已经把编译的字节码放在__pycache__文件夹中，这样以后再次运行的话，如果被调用的模块未发生改变，那就直接跳过编译这一步，直接去_pycache_文件夹中去运行相关的 *.pyc 文件，大大缩短了项目运行前的准备时间。

基本事实：

阅读代码的时间远大于写代码的时间
学习编程这类实践类的知识就是基本的就是试错
不谋全局着不足以谋一隅

scanpy源码浅析：或谈python面向对象结构

`_ 的含义`

`__ 的含义`

相关文章

网友评论

延伸阅读

深度阅读

栏目导航

热点阅读

单细胞分析之Scanpy

python单细胞分析

Single cell analysis

python学习收藏夹