6. Getting started with anndata

https://anndata.readthedocs.io/en/latest/tutorials/notebooks/getting-started.html

(更新日:2024-2-6 with ShortCake v1.9.1)

AnnData (Annotated Data) はScanpyを始め数多くのPython系シングルセル解析ツールで用いられているデータオブジェクト形式です。 AnnDataは行列形式のデータのために設計されています。すなわちn個のオブザベーションがあり、それぞれがd次元ベクトルとして表現され、各次元が変数または特徴に対応することを意味する。この行列の行と列は、インデックスが付けられます。

例えばscRNA-seqデータでは、各行が細胞に対応し、各列が遺伝子に対応する。さらに、各細胞と各遺伝子に対して、(1)各細胞のドナー情報、(2)各遺伝子のIDの代わりのシンボルのような付加的なメタデータがあるかもしれない。さらに、プロットに使用するカラーパレットのような構造化されていないメタデータがあるかもしれない。Pythonベースの他のデータ構造で以下の全てを満たすものは他に存在しない: - 疎な行列を扱える - 構造化されていないデータを扱える - 行と列それぞれのメタデータを扱える - user-friendly

[1]:
import numpy as np
import pandas as pd
import anndata as ad
from scipy.sparse import csr_matrix
print(ad.__version__)
0.8.0

まず、遺伝子発現量を格納する基本的なAnnDataオブジェクトを構築することから始めよう。

[2]:
counts = csr_matrix(np.random.poisson(1, size=(100, 2000)), dtype=np.float32)
adata = ad.AnnData(counts)
adata
[2]:
AnnData object with n_obs × n_vars = 100 × 2000

AnnDataがデータの要約を提供していることがわかる。

ここで与えた初期データは、adata.X を使ってスパース行列としてアクセスできます。

[3]:
adata.X
[3]:
<100x2000 sparse matrix of type '<class 'numpy.float32'>'
        with 126692 stored elements in Compressed Sparse Row format>

ここで、行列のobsvarそれぞれにインデックスを与えます。これには.obs_names.var_namesを使います。

[4]:
adata.obs_names = [f"Cell_{i:d}" for i in range(adata.n_obs)]
adata.var_names = [f"Gene_{i:d}" for i in range(adata.n_vars)]
print(adata.obs_names[:10])
Index(['Cell_0', 'Cell_1', 'Cell_2', 'Cell_3', 'Cell_4', 'Cell_5', 'Cell_6',
       'Cell_7', 'Cell_8', 'Cell_9'],
      dtype='object')

6.1. AnnDataのサブセット

これらのインデックス値はAnnDataをサブセットするために使用できる。 たとえば、特定の細胞クラスターや遺伝子モジュールの抽出に利用できる。 AnnData をサブセットするルールは、Pandas の DataFrame とよく似ている。 obs/var_names の値、ブーリアン・マスク、セル・インデックスの整数を使用することができる。

[5]:
adata[["Cell_1", "Cell_10"], ["Gene_5", "Gene_1900"]]
[5]:
View of AnnData object with n_obs × n_vars = 2 × 2

6.2. 整列されたメタデータの追加

6.2.1. Observation/variableレベル

オブジェクトのコアができたので、次はObservationとvariableにメタデータを追加しましょう。 adata.obsadata.varはPandas DataFramesであるため、これは簡単です。

[6]:
ct = np.random.choice(["B", "T", "Monocyte"], size=(adata.n_obs,))
adata.obs["cell_type"] = pd.Categorical(ct)  # Categoricals are preferred for efficiency
adata.obs
[6]:
cell_type
Cell_0 B
Cell_1 B
Cell_2 B
Cell_3 Monocyte
Cell_4 B
... ...
Cell_95 T
Cell_96 B
Cell_97 T
Cell_98 B
Cell_99 T

100 rows × 1 columns

[8]:
# 追加されたメタデータを確認
adata
[8]:
AnnData object with n_obs × n_vars = 100 × 2000
    obs: 'cell_type'

6.2.2. メタデータを使ってサブセット

[9]:
bdata = adata[adata.obs.cell_type == "B"]
bdata
[9]:
View of AnnData object with n_obs × n_vars = 26 × 2000
    obs: 'cell_type'

6.3. Observation/variableレベルの行列

また、UMAPエンベディングのように、多くの次元を持つメタデータをいずれかのレベルで持つかもしれない。 このようなメタデータのために、AnnDataは.obsm/.varm属性を持ちます。追加する異なる行列を識別するためにkeyを使用する。 .obsm/.varmはそれぞれ独立に異なる次元数を持つことができるが、.obsm行列は.n_obsと等しい長さでなければならず、.varm行列は.n_varsと等しい長さでなければならない。

保存したいデータのUMAP2次元エンベディングとして解釈できるランダム行列と、遺伝子レベルのいくつかのランダムメタデータから始めてみよう:

[10]:
adata.obsm["X_umap"] = np.random.normal(0, 1, size=(adata.n_obs, 2))
adata.varm["gene_stuff"] = np.random.normal(0, 1, size=(adata.n_vars, 5))
adata.obsm
[10]:
AxisArrays with keys: X_umap
[11]:
# adataの表示もアップデートされている
adata
[11]:
AnnData object with n_obs × n_vars = 100 × 2000
    obs: 'cell_type'
    obsm: 'X_umap'
    varm: 'gene_stuff'

.obsm/.varmに関するいくつかの注意点

  • 「配列のようなメタデータ」は、Pandas DataFrame、scipy sparse matrix、numpy dense arrayから得ることができます。

  • Scanpyを使用する場合、それらの値(列)は簡単にプロットされません。代わりに、.obsからの項目は、例えばUMAPプロットに簡単にプロットされます。

6.4. 構造化されていないメタデータ

AnnDataには.unsがあり、構造化されていないメタデータを入力することができる。 たとえばデータの分析に役立つ一般的な情報を含むリストや辞書のようなものが追加できる。

[13]:
adata.uns["random"] = [1, 2, 3]
adata.uns
[13]:
OverloadedDict, wrapping:
        OrderedDict([('random', [1, 2, 3])])
With overloaded keys:
        ['neighbors'].
[14]:
adata
[14]:
AnnData object with n_obs × n_vars = 100 × 2000
    obs: 'cell_type'
    uns: 'random'
    obsm: 'X_umap'
    varm: 'gene_stuff'

6.5. レイヤー

最後に、元のコア・データの異なる形式、たとえば正規化された発現量とされていない発現量があるかもしれない。 これらはAnnDataの異なる「レイヤー」に格納することができます。 例えば、元のデータをログ変換して違うレイヤーに保存してみましょう:

[15]:
adata.layers["log_transformed"] = np.log1p(adata.X)
adata
[15]:
AnnData object with n_obs × n_vars = 100 × 2000
    obs: 'cell_type'
    uns: 'random'
    obsm: 'X_umap'
    varm: 'gene_stuff'
    layers: 'log_transformed'

6.5.1. データフレームへの変換

レイヤーの1つをDataFrameとして返すこともできる:

[16]:
adata.to_df(layer="log_transformed")
[16]:
Gene_0 Gene_1 Gene_2 Gene_3 Gene_4 Gene_5 Gene_6 Gene_7 Gene_8 Gene_9 ... Gene_1990 Gene_1991 Gene_1992 Gene_1993 Gene_1994 Gene_1995 Gene_1996 Gene_1997 Gene_1998 Gene_1999
Cell_0 0.693147 0.693147 0.693147 0.693147 0.693147 0.693147 0.693147 0.693147 0.693147 0.693147 ... 1.098612 0.000000 0.000000 0.693147 0.693147 0.000000 0.000000 1.098612 0.000000 0.693147
Cell_1 1.098612 0.693147 0.000000 0.000000 0.693147 0.693147 0.693147 0.693147 0.000000 0.000000 ... 1.386294 0.693147 0.000000 0.000000 0.000000 0.000000 0.693147 1.098612 0.000000 0.693147
Cell_2 0.693147 0.693147 1.098612 0.000000 0.693147 0.693147 0.693147 0.693147 0.693147 0.000000 ... 0.000000 0.000000 0.000000 0.693147 0.000000 1.386294 0.693147 0.693147 0.693147 1.386294
Cell_3 0.000000 0.000000 0.693147 1.386294 0.693147 0.000000 0.693147 1.098612 0.693147 0.693147 ... 0.000000 0.000000 0.000000 0.000000 1.609438 0.000000 1.098612 1.386294 1.386294 0.693147
Cell_4 0.693147 0.693147 0.000000 0.693147 0.693147 0.693147 0.693147 0.000000 0.000000 0.693147 ... 1.098612 0.693147 1.098612 0.693147 1.386294 1.098612 1.098612 0.693147 1.386294 0.693147
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
Cell_95 0.000000 0.693147 0.693147 0.693147 0.693147 1.098612 0.000000 1.098612 1.386294 0.000000 ... 0.693147 0.693147 0.000000 0.693147 0.693147 0.693147 0.000000 0.693147 1.098612 0.000000
Cell_96 0.000000 0.000000 0.693147 0.693147 1.098612 0.693147 0.000000 0.000000 0.000000 0.693147 ... 1.098612 0.693147 1.609438 0.000000 0.693147 0.693147 1.098612 0.693147 0.693147 1.098612
Cell_97 0.693147 1.098612 0.693147 0.693147 0.000000 0.693147 0.693147 0.000000 0.000000 0.693147 ... 0.000000 0.693147 0.693147 0.000000 0.000000 1.098612 1.098612 1.098612 0.693147 0.693147
Cell_98 0.000000 0.000000 0.693147 0.000000 0.000000 0.000000 0.693147 1.098612 0.000000 0.693147 ... 0.000000 0.000000 1.098612 1.098612 0.000000 0.693147 1.098612 0.000000 0.693147 0.000000
Cell_99 0.000000 0.000000 1.098612 0.693147 0.000000 0.693147 1.386294 0.693147 0.000000 1.098612 ... 0.000000 0.000000 1.386294 0.000000 0.000000 1.386294 1.386294 1.098612 0.000000 0.000000

100 rows × 2000 columns

.obs_names/.var_namesがこのPandasオブジェクトの生成に使用されていることがわかる。

6.6. 結果をディスクに書き込む

AnnDataは、HDF5ベースの永続的な独自ファイル形式.h5adを提供する。 カテゴリ数が少ない文字列列がカテゴリカル形式でない場合、AnnDataは自動的にカテゴリカルに変換する。

[17]:
adata.write('my_results.h5ad', compression="gzip")
[18]:
!h5ls 'my_results.h5ad'
X                        Group
layers                   Group
obs                      Group
obsm                     Group
obsp                     Group
uns                      Group
var                      Group
varm                     Group
varp                     Group

6.7. まとめ

AnnDataはPythonによるシングルセル解析のスタンダードとなっているが、それには理由がある。使い方が簡単で、キーベースのストレージにより再現性の高い分析が可能になる。 シングルセル解析のための一般的なRベースのフォーマットへの変換も容易になってきている。

[19]:
import session_info
session_info.show()
[19]:
Click to view session information
-----
anndata             0.8.0
numpy               1.24.3
pandas              2.0.3
scipy               1.10.1
session_info        1.0.0
-----
Click to view modules imported as dependencies
anyio                       NA
asciitree                   NA
asttokens                   NA
attr                        23.1.0
attrs                       23.1.0
babel                       2.11.0
backcall                    0.2.0
bottleneck                  1.3.5
brotli                      1.0.9
certifi                     2024.02.02
charset_normalizer          2.0.4
cloudpickle                 3.0.0
colorama                    0.4.6
comm                        0.1.2
cython_runtime              NA
dask                        2023.5.0
dateutil                    2.8.2
debugpy                     1.6.7
decorator                   5.1.1
entrypoints                 0.4
executing                   0.8.3
fasteners                   0.19
fastjsonschema              NA
h5py                        3.8.0
idna                        3.4
importlib_metadata          NA
importlib_resources         NA
ipykernel                   6.25.0
ipython_genutils            0.2.0
jedi                        0.18.1
jinja2                      3.1.3
json5                       NA
jsonschema                  4.19.2
jsonschema_specifications   NA
jupyter_server              1.23.4
jupyterlab_server           2.25.1
markupsafe                  2.1.3
mkl                         2.4.0
mpl_toolkits                NA
natsort                     8.4.0
nbformat                    5.9.2
numcodecs                   0.12.1
numexpr                     2.8.4
packaging                   23.1
parso                       0.8.3
pexpect                     4.8.0
pickleshare                 0.7.5
pkg_resources               NA
platformdirs                3.10.0
prometheus_client           NA
prompt_toolkit              3.0.43
psutil                      5.9.0
ptyprocess                  0.7.0
pure_eval                   0.2.2
pyarrow                     15.0.0
pydev_ipython               NA
pydevconsole                NA
pydevd                      2.9.5
pydevd_file_utils           NA
pydevd_plugins              NA
pydevd_tracing              NA
pygments                    2.15.1
pytz                        2023.3.post1
referencing                 NA
requests                    2.31.0
rfc3339_validator           0.1.4
rfc3986_validator           0.1.1
rpds                        NA
send2trash                  NA
setuptools                  68.2.2
setuptools_scm              NA
six                         1.16.0
sniffio                     1.3.0
socks                       1.7.1
stack_data                  0.2.0
terminado                   0.17.1
tlz                         0.12.1
tomli                       2.0.1
toolz                       0.12.1
tornado                     6.3.3
traitlets                   5.7.1
typing_extensions           NA
urllib3                     1.26.18
wcwidth                     0.2.5
websocket                   0.58.0
yaml                        6.0.1
zarr                        2.16.1
zipp                        NA
zmq                         23.2.0
zope                        NA
-----
IPython             8.12.2
jupyter_client      7.4.9
jupyter_core        5.5.0
jupyterlab          3.6.3
notebook            6.5.4
-----
Python 3.8.18 (default, Sep 11 2023, 13:40:15) [GCC 11.2.0]
Linux-5.15.0-91-generic-x86_64-with-glibc2.17
-----
Session information updated at 2024-02-06 22:22
[ ]: