由于生信数据下载实在是太麻烦了!!网上的方法太多太冗杂,特写篇博客记录一下最近半年用到的方法。值得注意的是,以下部分代码不是原创(不过整理是原创的!),但是出处已经找不到了。此篇博客仅供学习交流作用。

1.单细胞数据获取

1.1借助包来下载GEO官网的数据

1.1.1方法一:利用GEOquery

        以GSE213252为例,可以用GEOquery来做到。

if(!requireNamespace('BiocManager',quietly = TRUE))
  + install.packages('BiocManager')
BiocManager::install('GEOquery')
library(GEOquery)
cancer <- getGEO('GSE213252')
cancer1 <- cancer[[1]]

1.1.2方法二:利用GEOquery

        以GSE148673为例,如果方法一行不通,可以试试运行下方代码。

suppressMessages(library(GEOquery))
library(Biobase)
options( 'download.file.method.GEOquery' = 'libcurl' ) 
gds <- getGEO("GSE148673" , destdir = "./",AnnotGPL = FALSE,getGPL = FALSE)
gds1 <- gds[[1]]

1.1.3方法三:更改国内镜像

  #the function used in this method require "a file"
  #"Your GSE may not be expression by array, or even not a GSE"
install.packages("devtools")
devtools::install_git("https://gitee.com/jmzeng/GEOmirror")
remotes::install_github("jmzeng1314/AnnoProbe")
library(GEOquery)
library(GEOmirror)
library(AnnoProbe)
eSet=AnnoProbe::geoChina('GSE131928')
eSet
eSet=eSet[[1]]
phenoDat <- pData(eSet)

1.2手动下载并导入数据

        顾名思义,先去GEO官网下载对应数据,然后再导入r或python里。以下根据不同的数据格式提供多种方法参考。

1.2.1方法一:features\barcodes\matrix

        无论是r还是python都要记得把这三个文件都改名,并放在同一个文件夹下。如下图所示。

【r版导入成Seurat对象】
library(Seurat)
setwd("/Users/mhuang/Downloads")
dir <- "./data/"
list.files(dir)
counts <- Read10X(data.dir = dir)
class(counts)
scRNA <- CreateSeuratObject(counts = counts)
scRNA
colnames(scRNA@meta.data$cell_types) #按理说是放在这里的
【python版导入成adata】

        python直接用scanpy包的read_10x_mtx函数导入即可。以我上方的文件夹截图为例,由于我是control_rep1文件夹里存放着这三个文件,因此路径只要写到control_rep1就好。

import scanty as sc
adata = sc.read_10x_mtx("/Users/mhuang/Downloads/control_rep1")

        如果python版导入的时候出现问题了,那可能是文件本身不是标准格式,是需要更改的。以我上方的文件夹截图为例,由于文件太大, 因此都是以gz结尾的压缩包。如果我导入失败,那么我首先应该在终端上解压文件(以matrix文件为例,解压方式如下图所示),然后检查它是否符合标准格式,若不符合,就更改。

gunzip -c matrix.mtx.gz | head

标准格式可以查看下方的参考文献:

seurat读取文件的格式 10x文件内容 mtx格式scanpy.read_10x_mtx scanpy读取10x格式_10x单细胞mtx文件中的基因列和细胞列互换了,如何换回来-CSDN博客

1.2.2方法二:csv

        方法二和方法三单纯是我记不住函数,所以全部都写出来一下!

xb=read.csv("/Users/mhuang/Downloads/counts.table.csv")
   #两个参数的设置和read.table函数一样

1.2.3方法三:txt\tsv

file=read.table("/Users/mhuang/Downloads/counts.table.tsv",fill = TRUE,sep="\t")
                                              #如果有缺失值就可以加上这个指标
colnames(file)=file[1,] #将第一行作为列名(或者read.table加上参数header=TRUE)
rownames(file)=file[,1] #将第一列作为行名(或者read.csv加上参数row.names=1)
file <- file[-1,] #删除第一行

1.2.4方法四:h5\hdf5

        这边仍然提供两种方法,如下所示。

  #1.
BiocManager::install("rhdf5")
library(rhdf5)
h5_file= H5Fopen("new.h5")
####new.h5文件内创建了一个组(group1_mat),组内又创建了df和matrix两个层级用以保存矩阵和数据框
h5dump(h5_file,load=FALSE)

  #2.
BiocManager::install("rhdf5")
library(rhdf5)
library(Seurat)
data_sample <- Read10X_h5("xxx.h5") 
data_seurat <- CreateSeuratObject(data_sample,project = "data_sample")

2.单细胞数据转化

2.1h5ad转为Seurat

2.1.1方法一:利用seurat-disk包来转化

        这个方法是最常见的教程,但是我基本上没运行成功过。。

remotes::install_github("mojaveazure/seurat-disk")
library(SeuratDisk)
Convert("xxxxx.h5ad", dest="h5seurat",
        assay = "RNA",
        overwrite=F)
seurat_object <- LoadH5Seurat("xxxxx.h5seurat")

2.1.2方法二:先转为h5和tsv文件,再导入r中

        这个方法百试百灵,应该还能运用到st数据(空转数据)中,不过我没试过,有需要可以试试。主要是分为两步:

第一步,把注释信息和表达矩阵在python中分离出来。

import h5py
import pandas as pd
import scanpy as sc
adata = sc.read("/file_root/PBMC_Batch_test.h5ad")
data_R=pd.DataFrame(data=adata.X,  #本来要加个todense()的,但是这边不知道为啥已经是密集矩阵了
                    index=adata.obs_names,
                    columns=adata.var_names)
data_R.to_hdf("/file_root/process_test_res.h5","data_R")
meta_data=pd.DataFrame(data=adata.obs)
meta_data.to_csv('/file_root/process_test_meta.tsv',sep="\t") 

第二步,导入进r并创造seurat对象。

data_R <- h5read("/file_root/process_test_res.h5","data_R")

mat <- data_R$block0_values
rownames(mat) <- data_R$axis0
colnames(mat) <- data_R$axis1
mat <- Matrix(mat, sparse = TRUE) #转为稀疏矩阵
meta_data <- read.table('/file_root/process_test_meta.tsv',sep="\t",header=T,row.names=1)

facs_test <- CreateSeuratObject(mat,project='Spatial',meta.data=meta_data)

2.2SingleCellExperiment对象转为h5ad

        主要是用python实现的转化,其中利用rpy2包来运行r代码。

import numpy as np
import pandas as pd
import scanpy as sc
import anndata as ad
from rpy2.robjects.packages import importr

importr("Seurat")
r('''
sce = readRDS('dataset_C(1).Rds')    
    ''')
obs = pd.DataFrame()
celltype = r('sce@colData@listData[["X"]]')
obs['cell.type']= pd.DataFrame(celltype)
X = r('sce@assays@data@listData[["counts"]]')
X = np.array(X).T
index = r('sce@colData@rownames')
adata = ad.AnnData(X, obs=obs, dtype='int32')
adata.obs.index = index
adata.write_h5ad('adata4.h5ad')

2.3 SingleCellExperiment对象转为Seurat

library(Seurat)
colnames(pbmc.sce) <- make.unique(colnames(pbmc.sce))#要确保列名唯一
pbmc.seurat <- as.Seurat(pbmc.sce)

2.4 Seurat转为h5ad

        其实是2.1.1节的函数反过来运行即可。

library(SeuratDisk)
library(Seurat)
SaveH5Seurat(seurat,filename = "dlpfc.h5Seurat")
Convert("dlpfc.h5Seurat",dest="h5ad")

Logo

腾讯云面向开发者汇聚海量精品云计算使用和开发经验,营造开放的云计算技术生态圈。

更多推荐