TCGA分析-数据下载-1

2023-11-01 20:20:00 浏览数 (2)

是https://cloud.tencent.com/developer/article/2353511 数据整理的上一步


title: "Untitled"

output: html_document

date: "2023-11-01"


R Markdown

代码语言:text复制
#实战代码有很多注意事项, 请不要不听课直接跑代码。
#数据下载
rm(list = ls())
library(GEOquery)
#加速的方法
#library(AnnoProbe)
#library(GEOquery)
#gse_number = "GSE56649"
#eSet=geoChina(gse_number)代替黑字前八行。
#先去网页确定是否是表达芯片数据,不是的话不能用本流程。
proj = "GSE218606"
eSet <- getGEO(proj, destdir = '.', getGPL = F)
代码语言:txt复制
## Found 1 file(s)
代码语言:txt复制
## GSE218606_series_matrix.txt.gz
代码语言:txt复制
## Using locally cached version: ./GSE218606_series_matrix.txt.gz
代码语言:text复制
class(eSet)
代码语言:txt复制
## [1] "list"
代码语言:text复制
length(eSet)
代码语言:txt复制
## [1] 1
代码语言:text复制
eSet = eSet[[1]]
#(1)提取表达矩阵exp
dat = rio::import("GSE218606_gene_count.txt.gz")
class(dat)
代码语言:txt复制
## [1] "data.frame"
代码语言:text复制
library(stringr)
k <- str_detect(colnames(dat), "^(NC|OMV2|gene_id|gene_name)") 
dat=dat[,k]
#或者也可以用dat=dat[,c(2:4,11:13)]
# 深坑一个
dat[97,3]
代码语言:txt复制
## [1] 19823
代码语言:text复制
as.character(dat[97,3]) #眼见不一定为实吧。庐山真面目
代码语言:txt复制
## [1] "19823"
代码语言:text复制
exp1 = dat
### 3.表达矩阵行名ID转换
#library(tinyarray)
#trans_exp_new是一个函数
#exp = trans_exp_new(exp)
#去重复的代码还可以是dat=distinct(dat,gene_name,.keep_all=T),.keep_all = T 可能是指定在删除重复项时是否保留所有信息。在某些情况下,当删除重复项时,可能会默认只保留第一行,而 .keep_all = T 可能指示保留所有重复行。但这取决于 distinct 函数的具体实现。
#dat=as.matrix(dat(,c[(2:4,11:13))])
#exp[1:4,1:4]
library(dplyr)
duplicates <- duplicated(exp1$`gene_name`)  
exp3 <- exp1[!duplicates, ]
rownames(exp3)=exp3$`gene_name`
exp3 <- exp3 %>% select(-`gene_name`)
exp3=exp3 %>% select(-`gene_id`)
#gdc下载的数据从此处开始衔接
### 4.基因过滤
##需要过滤一下那些在很多样本里表达量都为0或者表达量很低的基因。过滤标准不唯一。
#过滤之前基因数量:

# 3.基因过滤
##需要过滤一下那些在很多样本里表达量都为0或者表达量很低的基因。过滤标准不唯一。
#过滤之前基因数量:

#### 常用过滤标准1:
#仅去除在所有样本里表达量都为零的基因
exp33=as.matrix(exp3)
exp4 = exp33[rowSums(exp33)>0,]
nrow(exp4)
代码语言:txt复制
## [1] 27233
代码语言:text复制
#### 常用过滤标准2(推荐):
#仅保留在一半以上样本里表达的基因

exp5 = exp4[apply(exp4, 1, function(x) sum(x > 0) > 0.5*ncol(exp4)), ]
nrow(exp5)
代码语言:txt复制
## [1] 19333
代码语言:text复制
exp6 = exp5
#在R语言中,若要把fun应用到x的每一列,margin参数应该设置为1。 #1,函数会应用于矩阵的每一列(即,横向)。 #2,函数会应用于矩阵的每一行(即,纵向)。
#常用的过滤基因的标准

### 4.分组信息获取 一般使control在前 treat在后 要变成因子型 才具有顺序
#group_list=c("L","NC",each=4)
#\的意思是取消正则表达式 \d 在这里表示数字字符
g=str_split(colnames(exp6)," ",simplify = T)[,1]
group= g %>% str_remove("_\d")
table(group)
代码语言:txt复制
## group
##   NC OMV2 
##    3    3
代码语言:text复制
#在R语言中,使用factor(x, levels = c("NC", "OMV2"))会设定因子x的取值顺序为"NC"和"L"。也就是说,"NC"会排在"L"之前。
group=factor(group,levels=c("NC","OMV2"))
#str_split(colnames(exp3),"\(",simplify = T): str_split是一个字符串拆分函数,它将colnames(exp3)(即数据框或矩阵exp3的列名)按照"("(即左括号)进行拆分。simplify = T表示将结果简化为向量。
#[,2]: 这是一个子集操作符,用于从上一步的输出中提取第二个元素。
library(tinyarray)

#已经变成因子型变量,normal在前,tumor在后
table(group)
代码语言:txt复制
## group
##   NC OMV2 
##    3    3
代码语言:text复制
### 5.保存数据
save(exp6,group,proj,file = paste0(proj,".Rdata"))

学自生信技能树

0 人点赞