单细胞分析过程中的稀疏矩阵删减

2023-12-21 12:43:28 浏览数 (1)

引言

在单细胞转录组分析中,偶尔会出现电脑内存有限等情况,无法直接读取所有数据,这种时候可以考虑分析部分数据。

网上的教程提供了 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.gzfeatures.tsvbarcodes.tsvselected.tsv四个文件,其中selected.tsv为包含了所需细胞名的单列无表头行名的tsv文件。features 或 barcodes 文件如果没有解压,需要提前解压。

输出文件说明

输出文件被输出到results文件夹内。文件名分别为barcodes.tsv.gzfeatures.tsv.gzmatrix.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

结论

总而言之但是读进去了,但是也是真慢啊...

引用

  1. python 和 R 写出表达矩阵为稀疏矩阵 matrix.mtx.gz 的方法-CSDN 博客
  2. 「单细胞转录组系列」如何从稀疏矩阵中提取部分数据进行分析_单细胞稀疏矩阵-CSDN 博客

0 人点赞