之前我们对R语言中的数据类型和结构有了初步的认识:
数据类型:
数据结构:
在建立以上认知的基础上,我们可以尝试开始对生物医学研究中会接触到的数据进行处理,接下来将结合案例带大家一起完成这一部分的学习。
ps.本推文绝大多数代码来自生信技能树的小洁老师,少部分是自己手搓的
在工作目录下的文件可以直接读取, 例如, 若a.txt放在工作目录下, 则正确的读取代码如下:
read.table("a.txt")这个时候就有人要问了,主播主播,为什么我按照你的读取文件时出现了这样的的报错啊:
Error in file(file, “rt”) : cannot open the connection In addition:Warning message: In file(file, “rt”) : cannot open file ‘a.txt’: No suchfile or directory这个错误的意思是说,R语言无法找到a.txt这个文件,可能是因为文件不在当前的工作目录下,或者文件名拼写错误。解决这个问题的方法是确保a.txt文件确实存在于当前的工作目录中,并且文件名正确无误。
按道理来说,如果你先打开在你存放数据的文件夹里的.Rproj文件,工作目录是会默认设定好的,因此排除文件名错误之后可以用这个代码看看工作目录:
getwd()返回
[1] "C:/Users/Administrator/Desktop/生信技能树4月班/Day3/2026_04_09_R_day3/R_day3"如果工作目录不是你存放数据的文件夹,那么你可以使用setwd()函数来设置工作目录,例如:
setwd("D:/#your_folder_path#")#替换#your_folder_path#为你存放数据的文件夹路径如果工作目录下文件较多, 可以分类存放。将需读取的文件放入子目录, 只需在读取时添加子目录名称即可。
例如, 若a.txt放在工作目录下”raw”子目录下, 则正确的读取代码如下:
read.table("raw/a.txt") V1 V2 V3 1 KO-1 patient 45 2 KO-2 patient 31 3 KO-3 patient 47 4 KOC-1 normal 51 5 KOC-2 normal 40 6 KOC-3 normal 44绝对路径是指精确位置, 相对路径是指相对于当前工作目录的文件路径。
如果我们采用绝对位置,那如果我换个电脑、换个盘或者换个文件夹,代码就全部都要改,这个过程中但凡出了点纰漏可能整个就凉了,所以不推荐使用绝对路径,通过Rproj文件直接设定好工作目录它不香吗?
比如:
ex2 = read.csv("ex2.csv")不管我的Rproj文件放在哪里,只要ex2.csv文件和它在一个文件夹就永远能读取
刚才我们用到的read.table()和read.csv()等函数是R语言自带的,但是很多参数设置都非常恶心,那我们能不能绕开这个麻烦?
还真可以,data.table包中的fread函数读取速度快,智能, csv、txt、tsv等不同格式的数据可用相同代码读取,使用简单, 例如:
示例一:ex1.txt
library(data.table)ex1 = fread("ex1.txt")class(ex1) [1] "data.table""data.frame"OK,那现在为什么我一个变量有了两个数据结构呢?这是因为fread函数默认返回的是data.table数据结构,data.table是data.frame的一个增强版本,具有更高的性能和更多的功能,因此它同时具有data.frame和data.table的数据结构属性。
但是为了避免很多麻烦,我们可以 指定data.table = F, 产生的结果只是”data.frame”数据结构, 无需考虑data.table这种数据结构的特殊运算规则。 别问为什么,反正你直接默认加上就行了
ex1 = fread("ex1.txt",data.table = F)class(ex1) [1] "data.frame"示例二:ex2.csv
ex2 = fread("ex2.csv",data.table = F)#把第一列设置为行名,这个操作在Tidyverse介绍中已经讲过了ex2 = tibble::column_to_rownames(ex2,"V1")ex2 KO-1 KO-2 KO-3 KOC-1 KOC-2 KOC-3 DDX11L1 37 37 45 23 30 41 AL645608.8 5 4 1 7 12 8 LINC01786 8 10 5 4 11 5 TMEM240 10 17 10 16 20 10 MMP23B 6 10 8 7 10 16 CALML6 12 8 9 15 32 15 TMEM52 82 94 81 61 52 69 CFAP74 14 15 18 13 6 6 GABRD 76 71 59 37 22 27 PRKCZ-AS1 13 18 16 16 13 7示例三:ex3.txt
ex3 = fread("ex3.txt",data.table = F)ex3#SPOT_ID = V2 V3 1 ID locus.type SPOT_ID 2 probe1 NonCoding NR_120492 // RefSeq 3 probe2 NonCoding NR_027447 // RefSeq 4 probe3 Coding NM_130900 // RefSeq 5 probe4 Coding NM_001307930 // RefSeq 6 probe5 Coding NM_005259 // RefSeq 7 probe6 NonCoding NR_034148 // RefSeq 8 probe7 Coding NM_030967 // RefSeq 9 probe8 Coding NM_001289158 // RefSeq 10 probe9 Coding NM_145719 // RefSeq 11 probe10 NonCoding NR_110751 // RefSeq为啥原来的V1列没有被设置为行名呢?这是因为ex3.txt文件中前面有三行注释内容,fread函数默认会将这些注释内容作为数据的一部分进行读取,因此第一列的V1并没有被正确地识别为行名。
给刚才的的代码加个参数就可以了
ex3 = fread("ex3.txt",data.table = F,skip = 3)#用skip跳过注释内容所在的行ex3示例四:ex4.tsv
规则:行名不允许重复
ex4 = data.table::fread("ex4.tsv",data.table = F)ex4 = dplyr::distinct(ex4,symbol,.keep_all = T) #用distinct函数去除重复行,保留第一行,这个在Tidyverse教过了ex4 = tibble::column_to_rownames(ex4,"symbol") # 这个也讲过,但上次的例子是行变列名ex4 ENSEMBL sample1 sample2 sample3 sample4 TSPAN6 ENSG00000000003 5913 5122 9796 6589 DPM1 ENSG00000000419 4155 4751 782 7918 SCYL3 ENSG00000000457 793 3563 891 1373 C1orf112 ENSG00000000460 7767 1152 8879 4691 FGR ENSG00000000938 6725 4217 3208 8283 CFH ENSG00000000971 5824 3879 982 3545 MATR3 ENSG00000015479 8362 1978 6446 4760 POLR2J4 ENSG00000214783 1836 9201 2127 9084rio包用于简化数据导入和导出,号称”一行代码读写所有常见数据格式”。
这个函数很智能,可以根据扩展名自动选择读取方式,极大地简化读写操作。支持读取单个或多个工作簿的Excel文件。
示例一:ex1.txt
ex1.txt可以简单读取, 不需要额外指定参数。
ex1 = rio::import("ex1.txt")ex1 id group age 1 KO-1 patient 45 2 KO-2 patient 31 3 KO-3 patient 47 4 KOC-1 normal 51 5 KOC-2 normal 40 6 KOC-3 normal 44示例二:ex2.csv
无法设置行名, 需要用tibble::column_to_rownames另行设置。
ex2 = rio::import("ex2.csv")ex2 = tibble::column_to_rownames(ex2,"V1")ex2 KO-1 KO-2 KO-3 KOC-1 KOC-2 KOC-3 DDX11L1 37 37 45 23 30 41 AL645608.8 5 4 1 7 12 8 LINC01786 8 10 5 4 11 5 TMEM240 10 17 10 16 20 10 MMP23B 6 10 8 7 10 16 CALML6 12 8 9 15 32 15 TMEM52 82 94 81 61 52 69 CFAP74 14 15 18 13 6 6 GABRD 76 71 59 37 22 27 PRKCZ-AS1 13 18 16 16 13 7示例三:ex3.txt
ex3 = rio::import("ex3.txt",skip = 3)ex3没有设置注释的参数, 只能用skip=3跳过前3行。
ID locus.type SPOT_ID 1 probe1 NonCoding NR_120492 // RefSeq 2 probe2 NonCoding NR_027447 // RefSeq 3 probe3 Coding NM_130900 // RefSeq 4 probe4 Coding NM_001307930 // RefSeq 5 probe5 Coding NM_005259 // RefSeq 6 probe6 NonCoding NR_034148 // RefSeq 7 probe7 Coding NM_030967 // RefSeq 8 probe8 Coding NM_001289158 // RefSeq 9 probe9 Coding NM_145719 // RefSeq 10 probe10 NonCoding NR_110751 // RefSeq示例四:ex4.tsv
ex4 = rio::import("ex4.tsv")ex4 = dplyr::distinct(ex4,symbol,.keep_all = T)ex4 = tibble::column_to_rownames(ex4,"symbol")ex4 ENSEMBL sample1 sample2 sample3 sample4 TSPAN6 ENSG00000000003 5913 5122 9796 6589 DPM1 ENSG00000000419 4155 4751 782 7918 SCYL3 ENSG00000000457 793 3563 891 1373 C1orf112 ENSG00000000460 7767 1152 8879 4691 FGR ENSG00000000938 6725 4217 3208 8283 CFH ENSG00000000971 5824 3879 982 3545 MATR3 ENSG00000015479 8362 1978 6446 4760 POLR2J4 ENSG00000214783 1836 9201 2127 9084示例五:ex5.xlsx
用import函数读取该文件。代码如下:
ex5= rio::import("ex5.xlsx")ex5 id group age 1 KO-1 patient 45 2 KO-2 patient 31 3 KO-3 patient 47 4 KOC-1 normal 51 5 KOC-2 normal 40 6 KOC-3 normal 44示例六:ex6.xlsx
该文件有3个工作簿。用import_list函数读取该文件。代码如下:
ex6 = rio::import_list("ex6.xlsx")ex6$setosa Sepal.Length Sepal.Width Petal.Length Petal.Width Species 1 5.1 3.5 1.4 0.2 setosa 2 4.9 3.0 1.4 0.2 setosa$versicolor Sepal.Length Sepal.Width Petal.Length Petal.Width Species 1 7.0 3.2 4.7 1.4 versicolor 2 6.4 3.2 4.5 1.5 versicolor$virginica Sepal.Length Sepal.Width Petal.Length Petal.Width Species 1 6.3 3.3 6.0 2.5 virginica 2 5.8 2.7 5.1 1.9 virginica读入的结果, 数据结构为列表, 每个工作簿的内容是列表的一个元素。
导出的示例代码:
rio::export(ex1,file = "zz.csv")但是这并不意味着万事大吉了,在现实情况下可能会遇上一些神人文件,比如xxx.csv,结果用import()函数读出来各种出问题
结果仔细一看分隔符是分号而不是逗号,这个时候就需要指定分隔符了,可以在import函数里加一个format参数:
#ex7 = rio::import("xxx.csv",format =";")万事大吉
在这里我们主要涉及if语句、ifelse函数
其实它就是一个如果…那么…否则…的语句,纯判断
i = -1if (i < 0) print('down')if (i > 0) print('up')## [1] "down"这个代码意思是说,如果i小于0,就打印”down”,如果i大于0,就打印”up”,如果i等于0,则什么都不打印。i取-1的时候肯定是打印down
i = 1if (i > 0){ print('yes')} else { print("no")}## [1] "yes"这个代码是说如果i大于0就打印”yes”,否则打印”no”
i = 0if (i > 0){ print('yes')} else if (i < 0) { print('no')} else { print('0')}## [1] "0"这份代码是说如果i大于0就打印”yes”,如果i小于0就打印”no”,否则打印”0” 本质上是先判断i是否大于0,在这个基础上如果触发了第一个else,就进入第二个判断,即i是否小于0,纯纯的套娃
传统的if语句只支持单个逻辑值的判断。成批判断衍生了ifelse函数,可以接受多个逻辑值组成的向量。
i = 1ifelse(i > 0,"yes","no")## [1] "yes"i = c(-1,-2,1,2)ifelse(i > 0,"yes","no")## [1] "no" "no" "yes" "yes"这个代码是说,如果i大于0就打印”yes”,否则打印”no”, ifelse函数也可以进行套娃
在实际的数据操作中ifelse()+str_detect()可以用来根据样本信息把样本分组
library(stringr)samples = c("tumor1","tumor2","tumor3","normal1","normal2","normal3")k1 = str_detect(samples,"tumor");k1ifelse(k1,"tumor","normal")k2 = str_detect(samples,"normal");k2ifelse(k2,"tumor","normal")ifelse(k2,"normal","tumor")## [1] TRUE TRUE TRUE FALSE FALSE FALSE## [1] "tumor" "tumor" "tumor" "normal" "normal" "normal"## [1] FALSE FALSE FALSE TRUE TRUE TRUE## [1] "normal" "normal" "normal" "tumor" "tumor" "tumor"## [1] "tumor" "tumor" "tumor" "normal" "normal" "normal"在这个案例中,我们首先创建了一个包含样本名称的向量samples,其中包含了肿瘤样本和正常样本的名称。然后,我们使用str_detect()函数来检测样本名称中是否包含”tumor”或”normal”字符串,并将结果存储在逻辑向量k1和k2中,切记是和否的顺序不要搞反了
练习题
把下列样本分为Control和Vemurafenib两个组
library(stringr)title <- c( "A375 cells 24h Control rep1", "A375 cells 24h Control rep2", "A375 cells 24h Control rep3", "A375 cells 24h Vemurafenib rep1", "A375 cells 24h Vemurafenib rep2", "A375 cells 24h Vemurafenib rep3")title## [1] "A375 cells 24h Control rep1" "A375 cells 24h Control rep2" ## [3] "A375 cells 24h Control rep3" "A375 cells 24h Vemurafenib rep1"## [5] "A375 cells 24h Vemurafenib rep2" "A375 cells 24h Vemurafenib rep3"你可以先思考一下,看看能不能用ifelse函数结合str_detect()函数来完成这个任务。
如果是我我会这么做:
k3 <- str_detect(title,"Control");k3#先生成逻辑值group <- ifelse(k3,"Control","Vemurafenib") #用ifelse函数根据逻辑值分组,其实就相当于二分类了,如果要进行别的分析可以把它转化为因子table(group) #这里用table函数看看分组结果,看看有没有问题,我要验牌## [1] TRUE TRUE TRUE FALSE FALSE FALSE## group## Control Vemurafenib ## 3 3牌没有问题!但是我们再仔细思考一下,如果样本名是某个数据框的行名,如果我要标注样本所属的组别该怎么办呢?
还记不记得dplyr中的mutate函数?我们可以先把行名提取出来变成向量,然后用mutate函数结合ifelse函数来标注组别
dplyr包还有个优雅的功能,还记得我们刚才提到的多个条件的ifelse语句吗
ifelse(i>0,"yes",ifelse(i<0,"no","0"))## [1] "no" "no" "yes" "yes"dplyr给我们提供了一种类似于数学中分段函数的写法:
library(dplyr)case_when( i > 0 ~ "yes", i < 0 ~ "no", TRUE ~ "0" )## [1] "no" "no" "yes" "yes"这段代码的意思是说,如果i大于0就返回”yes”,如果i小于0就返回”no”,否则返回”0”,case_when函数可以让我们更清晰地表达多个条件的判断逻辑,避免了嵌套的ifelse函数带来的复杂性。
那么这个时候又有人问了,这个TRUE为什么这么突兀地出现了?我们尝试这样去理解:这个世界上永远为真的数据类型有且只有TRUE,TRUE永远为真。所以当前面两种情况都没满足的时候就把TRUE拿出来兜底了。
for循环的原理就是依次代数,即将变量x里的每个元素(代称i)带入到大括号里的代码中。
x <- c(5,6,0,3)for (i in x){ print(i)}## [1] 5## [1] 6## [1] 0## [1] 3这个东西我在第一堂课的笔记里就已经提到过
# 因为我全部装过了,所以最后只输出了包的名字pkglist <- c("data.table", "rio", "tibble", "ggplot2")for (i in pkglist){ print(i)if(!require(i,character.only = T))install.packages(i)}## [1] "data.table"## [1] "rio"## [1] "tibble"## [1] "ggplot2"apply 系列函数提供对向量、矩阵、数据框、列表的函数式操作,比for循环更简洁。
1. 创建示例数据
test <- as.matrix(iris[c(1,2,51,52),1:4])rownames(test) = NULLtest## Sepal.Length Sepal.Width Petal.Length Petal.Width## [1,] 5.1 3.5 1.4 0.2## [2,] 4.9 3.0 1.4 0.2## [3,] 7.0 3.2 4.7 1.4## [4,] 6.4 3.2 4.5 1.52. 标准用法apply()函数的三个参数分别是变量、维度和函数。维度可以取1和2两个值,其中1表示按行操作,2表示按列操作。
求test矩阵的每一列的中位数:
apply(test, 2, median)## Sepal.Length Sepal.Width Petal.Length Petal.Width ## 5.75 3.20 2.95 0.80求test矩阵的每一行的方差:
apply(test, 1, var)## [1] 4.750000 4.149167 5.622500 4.286667lapply 会对列表或向量的每个元素应用指定的函数,并返回一个列表。
1. 创建示例数据
test <- list(x = matrix( 36:33, 2), y = matrix( 32:35, 1), z = matrix( 30:27, 4))test## $x## [,1] [,2]## [1,] 36 34## [2,] 35 33## ## $y## [,1] [,2] [,3] [,4]## [1,] 32 33 34 35## ## $z## [,1]## [1,] 30## [2,] 29## [3,] 28## [4,] 272. 标准用法
对每个矩阵求平均值:
lapply(test,mean)## $x## [1] 34.5## ## $y## [1] 33.5## ## $z## [1] 28.5sapply 与 lapply 类似,但会尽量将结果简化为向量或矩阵。
对每个矩阵求平均值:
sapply(test,mean)## x y z ## 34.5 33.5 28.5计算每个元素的维度:
sapply(test,dim)## x y z## [1,] 2 1 4## [2,] 2 4 1如何用for循环批量读入工作目录下的所有数据文件(以.txt格式的文件为例)
这个代码是我自己手搓出来的,大家也可以试试别的思路
names <- list.files(pattern = "\\.txt$", ignore.case = TRUE)data_list <- list()for (i in names) { file_name <- names[i] data_list[[i]] <- read.table(file_name, header = TRUE, sep = "\t")}这一部分的内容相对于之前更加接近实际操作中我们需要去进行的处理:即导入事先编制好的大批量数据。我们不仅学会了经典语法中对数据进行读写和处理的方法,也学会了用其他开发者编写的R包来优化并赋能这个过程。在我的叙述中,我尽量用轻松的语气和类比的方法帮助自己和其他读者理解其中的一些概念和要点,其中或有谬误之处,还请大家包涵亦或是指正,也欢迎大家来交流。
ps.大家可能会发现,这里面的一些代码赋值是用=号,另一些则是用<-号。虽然绝大多数时候不讲究这么多,但是为了规范请尽量采用后者!--观点采自《R语言实战》
Thanks!
By 鼠鼠