引言
在单细胞转录组分析中,偶尔会出现电脑内存有限等情况,无法直接读取所有数据,这种时候可以考虑分析部分数据。
网上的教程提供了 python 和 R 两种代码1,2,但是实际操作中发现 R 代码并未提供正确的写出功能,所以本文以 python 作为示范。
文件层级
代码语言:tree复制.
├── main.py
├── data
│ ├── barcodes.tsv
│ ├── features.tsv
│ ├── matrix.mtx.gz
│ ├── selected.tsv
├── results
│ ├── barcodes.tsv.gz
│ ├── features.tsv.gz
│ └── matrix.mtx.gz
输入文件要求
输入文件应被放入data
文件夹内。需要matrix.mtx.gz
、features.tsv
、barcodes.tsv
和selected.tsv
四个文件,其中selected.tsv
为包含了所需细胞名的单列无表头行名的tsv
文件。features 或 barcodes 文件如果没有解压,需要提前解压。
输出文件说明
输出文件被输出到results
文件夹内。文件名分别为barcodes.tsv.gz
、features.tsv.gz
和matrix.mtx.gz
。输出文件可以被Seurat::Read10X
读入。
过程
代码语言:py复制from scipy.io import mmread
import pandas as pd
import numpy as np
# 读取表达矩阵
_index = pd.read_csv("./data/features.tsv", index_col=1, sep="t", header=None)
_index.index.name = None # 把索引列的列名去掉
_col = pd.read_csv("data/barcodes.tsv", sep="t", header=None)
_col.index.name = None # 把列名向量的名去掉
_data = mmread("data/matrix.mtx.gz").todense()
# 处理表达矩阵
## 挑选需求 col
_col.index
_selected = pd.read_csv("./data/selected.tsv", sep="t", header=None)
_selected.index.name = None # 把索引列的列名去掉
filtered_index = list(set(_col[0]).intersection(_selected[0]))
print(filtered_index.__len__())
# 加行名列名
rna_count = pd.DataFrame(
data=_data, index=_index.index, columns=_col.iloc[:, 0].tolist()
)
filtered_columns = rna_count[filtered_index]
rna_count = filtered_columns
print(rna_count.iloc[0:3, 0:2])
print("gene_ID_len : " str(rna_count.shape[0])) ### 获取表达矩阵基因数
print("cell_ID_len : " str(rna_count.shape[1])) ### 获取表达矩阵细胞数
# 重新写出 DataFrame 为 10X 格式的 sparse matrix 等相关文件
import os
import shutil
import gzip
import scipy
import time
fmt = "%Y-%m-%d %a %H:%M:%S"
Date = time.strftime(fmt, time.localtime(time.time()))
outdir = "results"
os.makedirs(outdir, exist_ok=True)
# Save matrix.mtx.gz
reAnno_count_sparse_mtx = scipy.sparse.coo_matrix(rna_count.values)
scipy.io.mmwrite(
os.path.join(outdir, "matrix.mtx"),
reAnno_count_sparse_mtx,
comment="Generate DateTime::" str(Date),
)
with open(os.path.join(outdir, "matrix.mtx"), "rb") as mtx_in:
with gzip.open(
os.path.join(outdir, "matrix.mtx") ".gz", "wb"
) as mtx_gz: # 创建一个读写文件'matrix.mtx.gz',用以将 matrix.mtx 拷贝过去
shutil.copyfileobj(mtx_in, mtx_gz)
os.remove(os.path.join(outdir, "matrix.mtx"))
# Save barcodes.tsv.gz
barcodesFile = pd.DataFrame(rna_count.columns)
barcodesFile.to_csv(
os.path.join(outdir, "barcodes.tsv.gz"), sep="t", header=False, index=False
)
# Save features.tsv.gz
featuresFile = _index
featuresFile.to_csv(
os.path.join(outdir, "features.tsv.gz"), sep="t", header=False, index=False
)
将文件保存为main.py
即可运行。
下面是用到的库。
代码语言:txt复制numpy==1.24.3
pandas==2.0.1
scipy==1.11.4
结论
总而言之但是读进去了,但是也是真慢啊...
引用
- python 和 R 写出表达矩阵为稀疏矩阵 matrix.mtx.gz 的方法-CSDN 博客
- 「单细胞转录组系列」如何从稀疏矩阵中提取部分数据进行分析_单细胞稀疏矩阵-CSDN 博客