注:本博客旨在分享个人学习心得,有不规范之处请多多包涵!


目录

  • Function 函数
  • data.frame 数据框
  • DNA序列读取与分析实例
  • 结束语


Function 函数

像其它的编程语言一样,R语言也有内置函数(如前面用到的c())和自定义函数。函数一般由三个重要的部分组成:输入参数,函数主体,返回参数。R语言的函数也是允许无输入参数或返回参数的的。以下的例子为在R语言中构造与调用函数:

#用R语言内置的function()函数来声明函数,并在括号内声明输入参数。可以用=为参数设置默认值
#getDouble为构造的函数的名称
getDouble <- function(num=1){ 
  result <- num * 2 #函数主体
  return(result) #返回参数用return()
}

myResult <- getDouble(num=20) #调用函数,myResult应为40

return()print()的区别就在与print()会在Console显示输出而return()不会。

data.frame 数据框

比起list,data.frame在计算生物学中R语言应用更为广泛。就像matrix一样,data.frame也是二维的数据结构,即由行列组成。然而,matrix的所有元素就像vector一样一定要是同一种数据类型,但data.frame只需保证每一列的数据类型相同。这就使得data.frame能够很方便地被用来读取类似于.csv.xlsx这样的可能有多种数据类型的表格文件。data.frame的构造请见下例:

id <- c("ENO2", "TDH3", "RPL39", "GAL4") #蛋白质的代号

protName <- c("enolase", "glyceraldehyde-3-phosphate dehydrogenase", "60S ribosomal protein L39", "regulatory protein Gal4") #蛋白质的全名

abundance <- c(24563, 22369, 16232, 32.3) #细胞中的含量,单位ppm

length <- c(437, 332, 51, 881) #蛋白质中氨基酸的数量

yeastProt <- data.frame(id, protName, abundance, length, stringsAsFactors = FALSE) 
#构造data.frame;每个输入vector的长度要相同
#vector的变量名会成为表头,vector的值会按列填入data.frame

构造完成的data.frame应如下图:

R语言 function 设置函数 r语言中function函数的作用_数据结构

data.frame的索引访问与matrix类似,具体见下例:

yeastProt[,1] #输出第1列的全部内容
yeastProt[2, ] #输出第2行的全部内容,包括表头
yeastProt[1, 2] #输出第1行第2列的元素,enolase
yeastProt[1:3, c(1, 3)] #输出第1-3行第1和3个元素,包括表头
yeastProt$length[1:2] #输出列标题为length的第1和第2个元素,437 332
yeastProt[1, "abundance"] #输出第1行列标题为abundance的元素,24563
yeastProt[yeastProt[, 4] > 500, 1] #输出yeastProt第4列大于500的行的第1列元素,GAL4

以下表格是vector, matrix, data.frame的互相转换:

vector转为

matrix转为

data.frame转为

vector

\

as.vector()

matrix

as.matrix()

\

as.matrix()

data.frame

as.data.frame()

as.data.frame()

\

在将其它数据类型转化成data.frame时,表头可能是不正确的。修改表头请见下例:

myMatrix <- matrix(data=1:12, ncol=3, nrow=4) #3x4的matrix
myDF <- as.data.frame(myMatrix) #转为data.frame
#转化完默认的列标题为V1, V2, V3... 默认行标题为1,2,3...
colnames(myDF) <- c("c1", "c2", "c3") #自定义列标题
rownames(myDF) <- c("r1", "r2", "r3", "r4") #自定义行标题

设置完成后,matrix被转变成了如下图的data.frame:

R语言 function 设置函数 r语言中function函数的作用_数据分析_02

DNA序列读取与分析实例

.fasta文件是常用的读取核酸或蛋白质序列产生的文件格式。文件中开头有一些说明信息,如来源物种和染色体编号;接下来就是以单个字母代表的核酸/蛋白质编码产生的序列。本次请从这个链接下载使用NC_045512.2版本COVID19的DNA测序数据,提取码: vfti
请看下例用R语言进行对这一序列的分析:

install.packages('seqinr', repos='http://cran.us.r-project.org') #安装第三方库
library("seqinr") #载入第三方库
covid <- read.fasta(file="covid19.fasta") #读取.fasta文件,注意路径最好写完整。读取完是list格式
covid_seq <- covid[[1]] #把完整序列单独存储为一个vector

DNA由4种base组成:A、T、C、G。A和T配对,C和G配对。根据查格夫法则(Chargaff’s law)可推导,DNA中A和T的数量相同,C和G的数量相同。因为C和G配对后由三条氢键链接,CG组合是不容易被拆散的,也就是DNA的高CG含量会减小它变异的可能性。
计算CG含量请见下例:

countTable <- table(covid_seq) #统计每个base出现的次数
cgCount <- countTable[["c"]]+countTable[["g"]] #C和G的数量
totalCount <- length(covid_seq) #序列的总长
cgContent <- cgCount/totalCount #CG在序列中的占比

另一个比较重要的概念是DNA word,即长度大于1个base的DNA组合。通过分析DNA word是否被过表达或低表达,可以在一定程度上推测该物种演化的轨迹。R语言 function 设置函数 r语言中function函数的作用_R语言 function 设置函数_03用来表示一个DNA word是否被过表达或低表达。对于一个长度为2的DNA word,计算公式如下:
R语言 function 设置函数 r语言中function函数的作用_R语言 function 设置函数_04
其中R语言 function 设置函数 r语言中function函数的作用_数据分析_05是DNA word的出现频率除以序列总长,R语言 function 设置函数 r语言中function函数的作用_生物信息学_06R语言 function 设置函数 r语言中function函数的作用_生物信息学_07是单个base的出现频率除以序列总长。计算某DNA word被如何表达请见下例:

baseCount <- count(covid_seq, 1) #用seqinr这个第三方库中的count()函数找长度为1的base的出现次数,和R内置的table()类似
wordCount <- count(covid_seq, 2) #找长度为2的DNA word的出现次数
totalLength <- length(covid_seq) #序列总长

#比如,我们要找TT这个DNA word的出现次数
rouTT <- (wordCount[["tt"]]/totalLength)/((baseCount[["t"]]/totalLength)*(baseCount[["t"]]/totalLength))

R语言 function 设置函数 r语言中function函数的作用_R语言 function 设置函数_03等于1时,这个DNA word为完全随机出现,与单个base出现概率的乘积相同。然而,当R语言 function 设置函数 r语言中function函数的作用_R语言 function 设置函数_03大于1时,它的出现概率大于随机出现的概率,说明它被过表达,所以可推测这个DNA word在这一物种的进化方面可能有某些有益作用。反之亦然。

结束语

下次会以分子生物学和生态学的数据为例,介绍R语言中的矩阵图、条形统计图、饼图、点阵图与线性回归、带状图的使用方法,敬请期待!有任何问题或想法欢迎留言和评论!