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>
ここで、行列のobs
とvar
それぞれにインデックスを与えます。これには.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.obs
とadata.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
[ ]: