1环境准备
搭建 Python 高效开发环境: Pycharm Anaconda
2安装 scanpy
代码语言:javascript复制pip install scanpy
3AnnData
1AnnData 介绍与结构
在做单细胞分析之前,先了解用于单细胞数据分析的一种数据结构--AnnData。
它一般作为 scanpy 单细胞分析包的数据存储格式。
主要由以下几部分构成:
功能 | 数据类型 | |
---|---|---|
adata.X | 矩阵数据 | numpy,scipy sparse,matrix |
adata.obs | 观察值数据 | pandas dataframe |
adata.var | 特征和高可变基因数据 | pandas dataframe |
adata.uns | 非结构化数据 | dict |
下面我们动手构建一个用于创建 AnnoData 的虚拟数据
代码语言:javascript复制import numpy as np
import pandas as pd
import anndata as ad
from string import ascii_uppercase
# 设置观测值数量
n_obs = 1000
# 生成观察时间
obs = pd.DataFrame()
obs['time'] = np.random.choice(['day 1', 'day 2', 'day 4', 'day 8'], n_obs)
# 设置特征名
var_names = [i*letter for i in range(1, 10) for letter in ascii_uppercase]
# 特征数量
n_vars = len(var_names)
# 特征注释数据框
var = pd.DataFrame(index=var_names)
# 生成数据矩阵
X = np.arange(n_obs*n_vars).reshape(n_obs, n_vars)
2AnnoData 初始化
代码语言:javascript复制# 初始化 AnnoData 对象
# AnnoData 对象默认使用数据类型为 `float32`, 可以更精确的存储数据
# 这里设置为整数,为了演示方便
adata = ad.AnnData(X, obs=obs, var=var, dtype='int32')
# 一般默认将变量或特征存储在数据框的行
# 查看数据
print(adata)
3AnnoData 切片特性
可以看到 AnnData 具有和 dataframe 或 Array 相似的长相,同样具备相似的特性,比如切片:
代码语言:javascript复制# 通过切片查看观测值和变量
print(adata.obs_names[:10].tolist())
print(adata.obs_names[-10:].tolist())
print(adata.var_names[:10].tolist())
print(adata.var_names[-10:].tolist())
# 查看矩阵
print(X)
4AnnoData 的 view 特性
AnnoData 可以实现与 numpy 中的 view 相似的功能。
换句话说就是,我们每次操作 AnnoData 时,并不是再新建一个 AnnoData 来存储数据,而是直接找到已经之前初始化好的 AnnoData 的内存地址,通过内存地址来直接改变 AnnoData 的值。这样做的好处是:
- 无需分配多余的内存
- 可以直接修改已经初始化后的 AnnoData 对象
view 可以使用 .copy() 来得到 AnnoData 对象。
代码语言:javascript复制# 查看 'A' 列的头三个元素
print(adata[:3, 'A'].X)
# 设置 'A' 列的头三个元素
adata[:3, 'A'].X = [0, 0, 0]
# 查看 'A' 列的头五个元素
print(adata[:5, 'A'].X)
其实我们在调用 .[] 时,AnnoData已经在内部实现了该操作,也就是说该 view 会成为保存数据的 AnnoData 对象。
但是,如果将 AnnoData 对象的 view 中的一部分赋值,该内容会复制一份并生成新的数据存储对象。
代码语言:javascript复制adata_subset = adata[:5, ['A', 'B']]
print(adata_subset)
adata_subset.obs['foo'] = range(5)
可以看到,这时赋值会直接将 AnnoData 对象复制一份。现在 adata_subset 会重新得到一块内存用于存储实际数据,而不再仅仅是对 adata 的内存地址引用。
5读取数据
代码语言:javascript复制import scanpy as sc
import pandas as pd
# 初始化数据
adata = sc.read(filename)
# 加入数据
anno = pd.read_csv(filename_sample_annotation)
# 加入样本分组信息
adata.obs['cell_groups'] = anno['cell_groups'] # categorical annotation of type pandas.Categorical
# 加入时间信息
adata.obs['time'] = anno['time'] # numerical annotation of type float
# 甚至可以直接赋值 dataframe
adata.obs = anno
6备份到本地
代码语言:javascript复制# 计算对象大小的函数
def print_size_in_MB(x):
print('{:.3} MB'.format(x.__sizeof__()/1e6))
# 查看 adata 对象大小
print_size_in_MB(adata)
# 查看是否备份
adata.isbacked
# 设置备份地址
adata.filename = './write/test.h5ad'
# 查看是否备份成功
adata.isbacked
可以看到,我们的 adata 对象已经备份成功,而且就在本地 './write/test.h5ad' 目录。
前边提到的 view 特性在这里同样适用,我们来看看 adata_subset 是否备份成功。
代码语言:javascript复制adata_subset.isbacked
adata_subset.filename = './write/adata_subset_test.h5ad'
adata_subset.isbacked
adata_subset 并没有被启用备份模式,重新设置备份模式。
需要注意的是:备份仅影响数据矩阵 X,所有注释信息都保留在内存中。如果想对全部数据的更改保存,则必须将导出到本地。
9保存数据
代码语言:javascript复制adata.write("./write/my_results.h5ad")
代码语言:javascript复制adata.write_csvs('./write/my_results_csvs', )
参考:https://anndata.readthedocs.io/en/latest/