美文网首页生物信息学生物信息学Bioinformatics
学习生物信息,我都应该掌握那些内容?

学习生物信息,我都应该掌握那些内容?

作者: superqun | 来源:发表于2019-03-08 10:32 被阅读117次
    目录

    持续更新ing.....
    最后更新时间:2019/3/8 10:28 AM

    1.linux命令

    1.1 linux的权限

    linux的权限

    权限:

    权限 意义 分数权制
    w 4
    r 2
    x 执行 1

    eg: 7=4+2+1 #代表w+r+x权限

    1.2 帮助文件

    man(详细讲解)和--help(使用语法简介)

    1.3 目录操作

    ls -r   #反序排列
    ls -sh  # 显示目录文件大小
    ls -d a*    #查看以a开头的目录或文件
    
    cd ~    #切换到当前用户目录
    cd -    #切换到上一次进入目录
    
    mkdir   #创建空目录
    mkdir -p    #创建多级目录
    rmdir -p    #删除空目录
    touch   #创建文件
    rm -r   #删除目录及目录下文件
    

    1.4 文件操作

    cp  #复制文件
    cp -r   #复制目录及目录中的文件
    mv  #移动/重命名 文件/目录:如果输入目录和输出目录一样则是改名字。
    

    1.5 文件编辑

    1.5.1 VI编辑器

    • 底行模式
      set nu:列出行号
      / charactor:查找字符
      wq:写入并退出
      q1:强制退出,不会保存修改
    cat     #查看文件内容,一次性输出
    cat 1.txt > 2.txt   #1.txt 文件内容取代2.txt文件
    cat 1.txt >> 2.txt  # 1.txt 文件插入2.txt文件末尾
    

    1.5.2 less/more

    less / more     #分页查看文件内容
    less -m     #显示百分比
    less -N     #显示行号
    /charactor  #查看模式下,向下查找字符。向上是 ?character
    b  #查看模式下下一页,下半页是d
    less -mN
    

    1.5.3 wc/head

    head / tail  #查看文件的头部/ 尾部
    wc  # 显示文件行数、字数
    wc -l <file>    #输出文件行数
    wc  -c <file>   #输出文件字节数
    wc -w <file>    #输出文件字节数
    wc -L <file>    #输出文件最长的一行
    ls |wc -l   #文件下有多少文件和目录
    
    ls |wc -l

    1.5.4 grep

    grep <str>  #文件查找,并输出匹配行
    grep -c <str>   #输出有几行匹配上了
    grep -i <str>   #不区分大小写
    grep -n <str>   #显示匹配行及其行数
    grep -v <str>   #不显示匹配行
    grep --color    #检索词加颜色,可以和-i -n 等一起用
    grep -n --color "the" 2.txt #结果如下
    
    cat 2.txt |grep -n --color "the"
    grep "t[ea]st" 2.txt    #匹配test 或者tast
    grep ^or$   #匹配行首或者行尾
    grep [^a]n  #匹配除an外含有n的行
    grep ^[TI]  #匹配以T或者I开头的行
    grep ^[^a-zA-Z] #不匹配字母开始的行
    
    在这里插入图片描述

    1.5.5 sort

    sort -n #依据数值大小排序
    sort -0 <file>  #排序结果存入制定文件
    sort -r #反序排列
    sort -k #以哪个区间排列
    sort -c #检查文件是否按照序列排序
    cat file1 file2 | sort | uniq > file3   #两文件并集
    cat file1 file2 | sort | uniq -d > file3    #两文件交集
    cat file1 file2 | sort | uniq -u > file3    #删除交集,留下各自补集
    

    1.5.6 sed

    sed -n  #列出经过sed处理过的一行
    sed -e  #直接在命令行模式进行sed动作编辑
    sed -i  #直接对原文编辑
    sed [参数] 'comand' 输入源文本  
    sed -n '1~2 p' a.txt    #a.txt第一行起始,2行为步长输出到屏幕(P)。奇数行。
    sed -n '1,+2 p' a.txt   #a.txt第一行起始,起始行后的连续+2行。即前三行。
    sed -n '/primary_transcript/,+4 p' genome.gff   #包含 primary_transcript 字符文本的一行及其后4行。
    sed '1,3 s/e/E/g' a.txt #对1-3行进行e到E的替换
    sed -e '1~2 d' a.txt    #删除奇数行
    sed -e '1,/^chr3/ d' a.txt  #定位:从1行开始到以chr3开头到行。操作:删除
    sed '^$ d' a.txt    #定位:空白行。操作:删除
    sed '2 a AT3GGO5780 GO:0000167' a.txt   #定位:第二行。操作:在第二行后追加一行。(在第二行之前是i)
    
    

    1.5.6.1 sed的定位

    • sed 第一行表示:1
    • 最后一行表示:$
    • 正则匹配:/REGEXP/
    • 步长处理:first~step
    • 选择选定行和后面几行:addr,+N

    1.5.6.2 sed的正则表达

    1.5.6.3 sed的操作命令**

    • a 在当前后添加一行
    • i 在当前行前添加一行
    • d 删除行
    • g 取缓存空间覆盖原有内容
    • p 打印行
    • ! 对选择之外对行执行sed
    • s 字符串替换

    1.5.7 cut

    cut -c 2-4 2.txt    #2.txt文件每行第2-4共三个字符截取
    cut -d '-' -f 2 2.txt   #2.txt文件以“-”为分割制表,每行第二列输出.-d '<charactor>'是指定分割符。默认是tab键
    

    在linux的shell中“tab”键表达方式是<kbd>control</kbd>+<kbd>V</kbd>+<kbd>I</kbd>

    1.5.8 文件操作-正则表达

    • 行首定位符:^ls | grep "^Green"grep "^RBM4" 2.txt
    • 行为定位符:$ls | grep "1$"
    • 限定符:
      *:前面字符出现0次或多次。>=0
      +:前面字符出现1次或多次。>=1
      ?:前面字符出现0次或一次。=0 OR =1
    • 字符集:
      [ ]:任意一个字符,数字和字母之间可用-连接表示范围。[abc][0-9][a-z][a-zA-Z]。(a | b | c) = [abc]。
      [^ ]: 字符非匹配。[^abc]不匹配a|b|c任意一个。

    正则表达练习网站:https://regexr.com/

    ━┉┉┉┉∞┉┉┉┉━━┉┉┉┉∞┉┉┉┉━━┉┉┉┉∞┉┉┉┉━━┉
    ●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞

    2. 测序结果的比对定量

    2.1 比对原理的介绍

    2.1.1 BWT算法学习

    需要进一步深入学习

    2.1.2 hisat2软件优势

    速度比较快,准确率高。

    2.2 bowtie2和hisat2的使用

    2.2.1索引构建

    hisat2-build Escherichia_coli.fa ./Escherichia_coli
    #hisat2对大肠杆菌基因组进行索引构建
    
    bowtie2-build Escherichia_coli.fa ./Escherichia_coli
    #bowtie2对大肠杆菌基因组进行索引构建
    

    2.2.2进行比对

    hisat2 -x Escherichia_coli -1 E_coli_1.test.clean.fq.gz -2 E_coli_2.test.clean.fq.gz  --un-conc-gz  ./E_coli_bowtie_unmap.fq.gz 2 > ./E_coli_bowtie_align.log   | samtools  sort -O BAM --threads 4 -o ./E_coli_bowtie_.bam
    # -x 指定基因组索引
    # -1 / -2  双端测序文件
    # --un-conc-gz 将没有匹配的序列存储到文件
    # >  将结果输出到log文件
    # -O  FORMAT 设置最终输出的文件格式,可以是bam,sam或者cram,默认为bam
    # -o  设置最后的文件名
    # --threads   设置线程数
    
    bowtie2 -x Escherichia_coli -1 E_coli_1.test.clean.fq.gz -2 E_coli_2.test.clean.fq.gz  --un-conc-gz  ./E_coli_hisat_unmap.fq.gz 2> ./E_coli_hisat_align.log   | samtools  sort -O BAM --threads 4 -o ./E_coli_hisat_.bam
    #注释同上
    

    2.3 计数

    2.3.1 why featurecount?

    1. featurecount对特征区域终止位置碱基处理更准确(与htseq相比)
    2. htseq对multmap处理更准确。a基因两reads匹配而b基因只有一个read匹配的话,则reads会匹配到a基因。而htseq则记录为multimap
    3. 更快。比htseq快22倍左右,内存占用为htseq的1/6

    2.3.2 count

    featureCounts -T 4 -F GTF -t exon -g gene_id -s 0 -Q 10 -C -B -p -a Escherichia_coli.exon.gtf -o ./feature_count.txt E_coli.bam
    # -T    线程数
    # -F    参考基因注释文件格式默认GTF
    # -t    指定gtf注释文件的feature类型。exon 即feature(count的对象)。默认将exon作为一个feature
    # -g    指定meta-feature类型。即从注释文件中提取Meta-features信息用于read count,默认是gene_id
    # -s    特异性建库
    # -Q    质量控制。控制比对reads最小质量值。默认0
    # -C    舍弃掉非正常匹配:paired read比对到不同染色体上或者相同染色体但不同链上
    # -B    paired reads都成功比对到参考基因组的reads才计数
    # -p    针对的是双端测序
    # -a    基因组注释文件GTF/GFF
    # -o    输出文件
    
    featureCounts -T 6 -p -t exon -g gene_id -a anno.gtf -o SRR_fCount.txt SRR.bam
    #官网的示例表达简写
    

    2.3.3 结果文件的提取

    cut -f 1,7 feature_count.txt |grep -v '^#' >featureCounts.txt
    # 提取第一列(Gene_id)和第七列(counts)数据。并去除空行。输出到文本文件
    

    下一步是差异基因分析。通过实验组和对照组的counts文件比对来分析DE:差异基因。
    然后是差异基因的选择。通过设置合适的差异基因筛选条件。padj<0.05,|log2FC|>=1(无生物重复)。
    筛选出合适的DE然后进行GO和KEGG分析。

    ━┉┉┉┉∞┉┉┉┉━━┉┉┉┉∞┉┉┉┉━━┉┉┉┉∞┉┉┉┉━━┉
    ●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞

    3 生物网站

    3.1 NCBI-BLAST

    blast访问地址

    在这里插入图片描述

    3.1.1 blast的使用

    如上图所示,blast分为blastn / blastp / blastx / tblastn / tblastx

    工具 查询序列 使用的数据库 适用范围
    blastn 核酸 核酸 核酸序列比对
    blastp 蛋白质 蛋白质 蛋白质序列与蛋白质数据库中的序列进行比较,可以寻找较远的关系
    blastx 核酸 蛋白质 将给定的核酸序列翻译成蛋白质与蛋白质数据库中的序列进行比对,对分析新序列和EST(表达序列标签)很有用
    tblastn 蛋白质 核酸 先将核酸数据库翻译为蛋白,在匹配。对于寻找数据库中序列没有标注的新编码区很有用
    tblastx 核酸 核酸 先将检索序列和核酸数据库全翻译为蛋白再匹配。特殊情况使用

    3.1.2 常用子数据库

    • NR数据库:非冗余蛋白库,所有 GenBank+EMBL+DDBJ+PDB中的非冗余蛋白序列,对于所有已知的或可能的编码序列, NR记录中都给出了相应的氨基酸序列以及专门蛋白 数据库中的序列号NR库相当于一个以核酸序列为基础的交叉索引,将核酸数据和蛋白数 据联系起来NT是NR的子集。(引自ppt)
      -** RefSeq数据库**:参考序列数据库,包含RefSeq_genomic, RefSeq_protein和RefSeq transpans,即具有生物意义上的非冗余基因,转录本和蛋白质序 列,是经过NCBI和其他组织校正的数据库,包括了官方的基因符号和可选的符号。(引自ppt)

    3.2 Ensembl数据库

    EMBL-EBI数据库包括Ensembl数据库(脊椎动物数据库)和Ensembl Genome(非脊椎动物全基因组数据库)数据库。

    Ensembl数据库中人类基因组和注释文件下载

    进入Ensembl首页点击上方导航栏的download进入下载界面,可以选择两种方式下载。可以点击通过FTP下载。进入/pub/目录。选择最新的数据文件。例如截止2019/3目前最新的是release-95/pub/release-95/下多种基因序列文件。基因组选择fasta/。注释文件可以选择gtf/。若要获取基因组文件,依次打开/pub/release-95/fasta/homo_sapiens/dna/选择相应文件下载。如果需要RNA-seq比对可以选择Homo_sapiens.GRCh38.dna_rm.toplevel.fa.gz下载。若要获取注释文件,依次打开/pub/release-95/gtf/homo_sapiens/需要RNA-seq回帖参考,可以选择Homo_sapiens.GRCh38.95.gtf.gz下载。
    </br>
    下载连接直达:
    人类基因组:ftp://ftp.ensembl.org/pub/release-95/fasta/homo_sapiens/dna/
    人类基因组注释文件:ftp://ftp.ensembl.org/pub/release-95/gtf/homo_sapiens/

    在这里插入图片描述 在这里插入图片描述

    3.3 UCSC

    随着众多物种基因组测序的完成,仅仅以纯文本的方式存储 和展示基因组DNA字符将无法满足对测序数据的研究,因此, UCSC应用而生。UCSC能够满足在任何尺度上快速地查询和显示 基因组的内容,以及对基因组序列进行注释,注释内容可以在一个 窗口中显示所有与某一基因组区域相关的信息:定位和序列信息、 已知基因和预测基因、表型和文献支持、mRNA和EST、调控 (CpG岛)、比较基因组信息、SNP、基因组重复元件等。(引自ppt)

    3.4 pfam-蛋白质结构预测

    需要进一步深入学习

    ━┉┉┉┉∞┉┉┉┉━━┉┉┉┉∞┉┉┉┉━━┉┉┉┉∞┉┉┉┉━━┉
    ●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞ ∞●∞

    4. R语言基础

    4.1 常用函数

    head()  #查看文件开头前十行.查看尾端信息是tail
    summary()   #统计总结。包括最小值,1/4位值,中位数,平均数,3/4位值,最大值
    dim()   #显示纬度。eg.输出结果是 150  5 代表有150行5列。
    model() #显示对象类型
    max()   #显示最大值,最小值是 min()
    sort()  #筛选排序。
    sort(iris$Sepal.Length,decreasing = T)  #iris数据的Sepal.Length列按照大到小排序(从小到大是decreasing = T)
    which.max() #返回最大值大索引下标。最小值是which.min()
    iris$Sepal.Length[which.max(iris$Sepal.Length)] #等价于max(iris$Sepal.Length)
    median()    #输出中位数
    length()    #显示对象中元素个数。length(iris)=5;length(iris$Sepal.Length)=150
    var()   #求数据的方差
    sum()   #对数据求和
    t()     #转换表格的坐标轴直接转换。
    
    
    

    4.2 数据类型

    R中的对象(数据集合):向量、矩阵、数据框、数组、列表。
    对象的五种基本类型:字符型(character)、数值型(numeric,包括小数 )、整型(integer)、复数型(complex)、逻辑型(logical,TRUE/FALSE)
    对象的属性:名字(names,dimnames),维度(dimensions, 包括矩阵等),类别(class,包括数字、整数等), 长度(length) attributes() 查看对象属性,没有属性返回NULL

    4.2.1 向量

    向量的创建

    a <- c(1,2,3,4) # 为a赋值numeric 数值型
    b <- c(1:4) #为b赋值integer数值型
    c <- c("hello","world") #赋值charac字符型
    d <- c(TURE,FALSE)  #赋值为logical逻辑型
    

    向量元素的访问

    a <- c(11:20)
    #生成以一个整形向量,11~20共十个元素
    a
    #列出a的所有元素
    a[1]
    #访问a中第一个元素即 11
    a[c(1,3,5,7)]
    #访问a中第1,3,5,7个元素,注意是 a[c(1,3,5,7)] ,而不是 a[1,3,5,7] 这是矩阵的访问格式
    a[3:8]
    #访问a中第3~8个元素
    
    在这里插入图片描述

    4.2.3 矩阵

    • 矩阵只能包含一种数据类型,如数值型。
    • 有二维性

    4.2.3.1 矩阵的创建

    x <- matrix(data = c(11:30) , nrow = 5 , ncol = 4 , byrow = FALSE)
    #matix  创建向量
    #data=c(11:30)  数据输入
    #nrow = 5,ncol = 4  有5行4列
    #byrow = FALSE  依据列填充,依据行填充是 byrow = TRUE
    
    y <- c(11:20)   #向量型y
    y   #查看向量型y内容
    dim(y)  #查看y的纬度信息
    dim(y) <- c(2,5)    #重新定义纬度信息2行5列
    dim(y)  #显示变量y的纬度信息
    y   #查看矩阵y的信息
    
    在这里插入图片描述
    z1 <- c(11:13)
    z2 <- c(21:23)
    z <- cbind(z1,z2)   #通过列对齐,对向量z1,z2合并.行对齐是rbind
    dim(z)
    
    在这里插入图片描述

    4.2.3.2 矩阵访问

    x <- matrix(data = c(11:30), nrow = 5, ncol = 4, byrow = FALSE)
    #创建矩阵
    x[3,3]
    #访问第3行第3列元素
    x[3,]
    #访问第3行所有元素
    x[,3]
    #访问第3列所有元素
    x[3]
    #访问第3个元素,顺序与矩阵生成byrow参数有关
    x[1:2,]
    #显示前两行元素
    x[,3:4]
    #显示第3到4列元素
    #------------------------------------------------------------------#
    #对于有标题到矩阵可以直接通过标题访问
    y <- matrix(data = c(11:30),nrow = 5,ncol = 4,byrow = FALSE,dimnames = list(c('r1','r2','r3','r4','r5'),c('c1','c2','c3','c4'))) 
    # 创建矩阵的过程中重命名行索引名称、列索引名称.如果不设置行名,可以填写NULL
    y['r3','c3']
    #访问第r3行r3列元素
    

    4.2.4 数据框(生信常用)

    • 可以存放多种数据的矩阵

    4.2.4.1 数据框的创建

    transcript_id <- c('TUCP_1','TUCP_2','lncRNA_1','lncRNA_2')
    # 创建字符型向量
    treat <- c(1,2,3,4)
    #创建数值型向量
    control <- c(5,6,7,8)
    # 创建向量
    chr <- c('chr1','chr2','chr3','chr4') 
    # 创建向量chr
    testdata <- data.frame(transcript_id,teat,ctrol,chr)
    # 创建数据框testdata
    dim(testdata) 
    # 查看数据框testdata的维度信息
    
    
    #--------------------#
    #也可以为数据库添加题头
    tdaddrowname <- data.frame(transcript_id,treat,control,chr,row.names = c('r1','r2','r3','r4'))
    # 为创建的数据框tdaddrowname更换行索引
    
    在这里插入图片描述

    4.2.4.1 数据框的访问

    taddrowname
    #查看数据框的内容
    tdaddrowname[,2:3] 
    #显示数据框中的第2列、第3列数据
    tdaddrowname[c('transcript_id','chr')] 
    #选取多列
    #矩阵访问的方式同样使用于数据框的访问
    tdaddrowname$transcript_id
    #如果数据框有列名,可以是直接使用$直接调用该列,直接调用该列
    
    
    #------------------------------------#
    inform <- names(tdaddrowname) %in% c('trancript_id','chr')
    #names()获取对象的或者设置对象的名称
    #%in%相当于match()函数的一个缩写。用来判断一个数组或矩阵是否包含在另一个数组或矩阵里。
    #A %in% B,若A中元素有和B中相匹配,则返回逻辑型 TURE,不匹配则返回逻辑型 FALSE
    newtd <- tdaddtowname[!inform]
    #反向选择非匹配的列
    newtd
    #查看
    

    4.2.5 数组

    • 数组通常是三纬的

    数组的创建

    x1 <- array(data = c(11:30))
    #创建一纬数组
    x2 <- array(data = c(11:30),dim = c(4,5))
    #创建二维数组
    length(x2)
    #查看元素个数
    dim(x2)
    #查看纬度信息
    x3 <- array(data = c(11:37),dim = (3,3,3))
    #创建三维数组
    x3 <- array(data = c(11:37),dim = c(3,3,3),dimname = list(c('r1','r2','r3'),c('c1','c2','c3'),c('m1','m2','m3')))
    #创建有索引标题的三维数组
    length(x3)
    dim(x3)
    x3
    #查看 x3 的属性和元素
    
    在这里插入图片描述

    数组的访问

    x1
    #一纬数组的访问和向量访问一样
    x2[2,3]
    x2[1,]
    #二维数组的访问和矩阵的访问方法一致
    x3['r1','c2','m3']
    #依据索引访问元素
    x3[1,2,3]
    #依据元素的纬度信息访问元素
    x3[,,3]
    #依据纬度信息,访问第三面的所有元素
    x3[,'c2',]
    #依据索引,访问所有面第二列元素
    
    在这里插入图片描述

    4.2.6 列表

    • 可以包含各种不同类型对象
    • 列表(list)内容中的代码引自ppt
    xlist <- list(c(1:5),matrix(data = c(1:10),nrow = 2,ncol = 5,byrow = FALSE),123,TRUE) 
    # 定义新列表
    length(xlist) 
    # 查看列表中的元素个数
    xlist
    # 查看列表xlist的全部信息
    dim(xlist) 
    # 列表没有维度信息
    dim(xlist) <- c(2,2)
    dim(xlist)
    class(xlist) 
    # 一旦给列表添加了维度信息,等价于转换成了矩阵格式
    
    ############# list访问 ############
    xlist <- list(c(1:5),matrix(data = c(1:10),nrow = 2,ncol = 5,byrow = FALSE),123,TRUE) 
    # 定义新列表
    xlist[[1]] 
    # 访问第一个元素
    class(xlist[[1]]) 
    # 属于向量中的范畴,可以在此基础上进行“向量访问”
    xlist[[1]][2:4] 
    # 访问列表第一个元素(向量类型)的第2至5的元素
    xlist[[2]][1,3] 
    # 访问列表中矩阵元素的第一行、第三列的元素
    xlist[[2]][2,] 
    # 访问列表中矩阵元素的第二行所有数据信息
    xlist[[c(2,4)]] 
    # 列表第二个元素的第4个元素
    

    4.3 数据处理

    4.3.1 数据导入

    • 常用read.table()函数
    xlsdata1 <- read.table(file = "FPKM.xls",header = TRUE, sep = '\t', stringsAsFactors = FALSE,row.names = 1)
    # read.table 读取文件
    # file="FPKM.xls"   要读取的文件
    ## header=TRUE  读取表头
    ## sep='\t' 设置列间隔是tab
    # stringsAsFactors = FALSE  #不作为因子?
    ## row.names = 1    设置行标题为第一列
    

    4.3.2 数据导出

    write.table(x = iris,file = 'test_iris1.xls',append = FALSE,sep = '\t',row.names = FALSE)
    # write.table() 数据导出
    # x=iris    数据集
    # file = 'test_iris1.xls'   导出文件名
    # append = FALSE    不以追加模式导出
    # sep='\t'  分隔符是TAB
    # row.name=FALSE    不导出行名
    

    4.3.3 缺省值处理

    • 缺省值会影响数据的处理
    • 对于NA值是错误、遗失的数据,NULL是不存在的,不占用空间,不会列于计算之中
    • 当向量数据中存在NA值,sum()结果也是NA值

    先定位再处理

    x <- c(1,3,NA,7)
    # 赋值向量
    is.na(x)
    ### 查询x中元素是否有NA值,如果是则在该元素纬度坐标中返回TURE,不是返回FALSE
    table(is.na(x))
    ### table()可以对元素进行统计并计数。
    sum(x,na.rm = TRUE)
    ## 在统计向量时候通过 na.rm = TRUE 来移除 NA 值,很多函数都有这个选项
    x[which(is.na(x))] <- 0
    ### which()函数返回符合条件的响应位置,which(is.na(x))会返回TURE对位置
    # which(b[,3]==1.7) 查询b的第三列中数值是1.7的位置
    # 把NA位置选出然后用0代替
    # 还可以使用 x[is.na(x)] <- 0 来代替NA
    
    在这里插入图片描述
    在这里插入图片描述
    其他方法处理NA值
    test <- data.frame(x=c(1,2,3,4,NA),y=c(6,7,NA,8,9))
    # 构建数据框
    
    which(is.na(test),arr.ind = T) 
    ### which(is.na(test)) 可以匹配NA值,arr.ind =T 可以返回相应返回行列坐标,适用于多维数据。如果没有arr.ind = T 则返回一个序列数。
    test[which(is.na(test),arr.ind = T)] <- 0 
    #结合which进行多维数据类型的缺失替代
    
    na.omit(test)
    ## 可以直接删除NA值所在的一行数据
    

    4.3.4 数据拆分与合并

    stack()

    data('PlantGrowth') #加载数据
    fomula(PlantGrowth) #查看格式
    
    pg <- unstack(PlantGrowth)
    # 将数据去堆栈排列
    stack(pg)
    # 将多列多数据堆栈排列
    stack(pg,select = ctrl)
    # 选择名为 ctrl 的特定向量
    stack(pg,select = -ctrl)
    # 选择并删除向量
    ## stack只使用于向量数据类型?
    

    cbind / rbind()

    x <- matrix(data = c(1:10),nrow = 5, ncol = 2)
    y <- matrix(data = c(21:30),nrow = 2,ncol = 5) 
    
    xyr <- rbind(x,y)
    # 将 x,y 依据行方向结合在一起
    xyc <- cbind(x,y)
    # 将x,y 依据列方向结合在一起
    ## cbind & rbind 结合的必须是相匹配的 例如cbind 需要列一致
    

    merge()

    ID <- c(1,2,3,4)
    name <- c('li','zhang','zhou','wu')
    score <- c(23,45,67,98)
    st1 <- data.frame(ID,name)
    st2<-data.frame(ID,score)
    
    st <- merge(st1, st2, by = 'ID')
    # 按照ID的对应关系将st1 st2的信息合并
    merge(ink1,ink2,by="id",all=T)  
    #所有数据列都放进来,空缺的补值为NA
    merge(ink1,ink2,by="id",all=F)  
    #默认,只取两者的共有的部分
    

    subset()

    filedata <- read.table(file = "FPKM.xls",header = TRUE,sep = '\t',stringsAsFactors = FALSE)
    class(filedata) 
    dim(filedata)
    summary(filedata)
    
    newset <- subset(x = filedata , subset = (treat2 > 100),select = treat1:control3)
    #subset 依据要求返回向量、矩阵、数据框的子集
    #x=filedata 数据输入
    #subset = (treat2 > 100)    选取treat2一列大于100的元素
    #select = treat1:control3   界定选取子集的 元素/范围/哪几列
    

    melt

    • reshape2包
    names <- c('A','B','C')
    yuwen <- c(87,67,98)
    shuxue <- c(78,69,88)
    yingyu <- c(56,64,77)
    xingbie <- c('F','M','F')
    chengjidan <- data.frame(names,yuwen,shuxue,yingyu,xingbie)
    chengjidan
    
    newtable <- melt(data = chengjidan,id.vars = c('names','xingbie'),variable.name = 'kemu',value.name = 'fenshu')
    

    4.3.5 数据的匹配

    R语言中有多种函数支持正则表达匹配,其中比较常用的是grep、gsub、grepl

    #完整的grep参数
    grep(pattern, x, ignore.case = FALSE, perl = FALSE, value = FALSE, fixed = FALSE, useBytes = FALSE, invert = FALSE)
    #grep() 查找特定格式元素。查找,存在参数value,返回结果是匹配项的下标
    #pattern    正则表达式
    #x  字符向量
    #ignore.case=FALSE  默认选项式FALSE表示区分大小写,不区分式TURE
    #perl=FALSE 默认false,表示兼容perl正则表达
    #fixed=FALSE    如果为TRUE,pattern是要匹配的字符串。覆盖所有冲突的参数
    #useBytes=FALSE 默认选项FALSE,TURE情况下是逐字节匹配而不是逐字符匹配
    #invert = FALSE 如果是TURE则不返回元素索引或值   
    
    函数 作用
    grep() 查找,存在参数value,返回结果是匹配项的下标
    grepl() 查找,返回值为true
    sub() 只对查找到的第一个内容进行替换。(同下)
    gsub() 对查找到的所有内容进行替换,返回替换后的text;否则 直接返回text

    表格引自ppt

    string1 <- 'qwewertyuioplkjhjnmbvc,fdsaa,email:chongqing@126.comsdderf'
    string2 <- '112351426325637,tel:0643-3445676,tel:0654-6778877,7655263529987352414099827tel:72625637784383726'
    
    grepl("[a-z]+@[0-9]+.c(n|om)",string1)
    #匹配 ‘字母@数字.com’或者'字母@数字.cn'格式的元素,有的话返回TURE
    grepl("[0-9]+-[0-9]+",string2)
    #匹配 ‘数字+-+数字’格式的元素
    

    5. R语言绘图

    R语言的四套图形系统:base、grid、lattice、ggplot2。随着ggplot2的发展越来越成为主流。

    5.1 绘图参数表示方式

    5.1.1利用元素属性直接更改

    常用?元素属性可以直接通过元素参数设置

    元素参数 意义
    font 字体
    lty 线类型
    lwd 线宽度
    pch 点类型
    xlab 横坐标
    ylab 纵坐标
    xlim 横坐标范围
    ylim 纵坐标范围

    5.1.2 属性.元素

    例如部分属性如下

    属性 意义
    col 颜色
    cex 缩放倍数
    font 字体

    部分元素如下

    元素 意义
    main 标题
    sub 副标题
    axis 坐标刻度
    lab 坐标轴标签

    因此。col.main 即 标题的颜色;cex.main 标题的缩放倍数;font.sub 副标题的字体样式。等等

    5.1.3 自定义属性需要用到函数修改

    一些特异性属性需要用到函数来修改
    legend(location, title, legend, ...)来修改图例的属性

    5.2 R语言绘图常用软件包概述

    5.2.1 画图软件

    画图 软件包
    热图 pheatap
    韦恩图 VeinnDiagram
    PCA图 ggplot
    箱线图/小提琴图/散点图/直方图/密度分布图/柱状图/气泡图/火山图 ggplot2

    5.3 ggplot2概述

    ggplot2特点

    • 图层设计
    • 统计变换
    • 扩展包丰富

    5.3.1 ggplot2 图形部件

    • 数据(data)和映射(mapping)
      反应的是数据/元素和图形对应关系
    • 标度(scale)
      标度负责控制映射后图形属性的显示方式 ,具体形式上来看是图例和坐标刻度(引自ppt)
    • 几何对象(Geometric)
      图像中的元素,控制绘制的图形类型
      在这里插入图片描述
    • 统计变换
      对原始数据进行某种计算,例如对二元散点图加上一 条回归线
    • 坐标系统
      坐标系统控制坐标轴并影响所有图形元素,坐标轴可 以进行变换以满足不同的需要
      可以通过变换坐标系统,绘制饼图、极坐标图
    • 图层系统
      数据、映射、几何对象、统计变换等构成一个图层
      图层可以允许用户一步步的构建图形,方便单独对图层进行修改
    • 分面
      控制分组绘图的方法和排列形式,将数据按某种方式分组,然后分别绘图

    5.3.2 ggplot2作图思路

    1. 数据映射
    2. 标度,图例
    3. 统计变换
    4. 图层,分面
    5. 美化

    5.3.3 图形元素

    在这里插入图片描述

    以上图为例,有以下元素

    • axis.title 坐标轴标题
    • axis.line 坐标轴线
    • axis.text 坐标轴标示
    • axis.text.bar 坐标轴bar
    • panel 面板
    • panel.grid 面板网格
    • panel.background 面板背景
    • legend 图例
    • legend.title 图例标题
    • legend.text 图例内容
    • legend.key

    5.4 R语言绘制柱状图

    library(ggplot2)
    setwd("~/SUPER QUN/3.2.2")  #设置工作目录
    fp <- read.table("FPKM.xls",header = T ,row.name = 1)   #读取目录
    head(fp)
    
    bar <- ggplot(fp) + geom_bar(aes(x = category, fill = strand),posirion = "fill")
    #ggplot(fp) 指定作图文件为fp
    #geom_bar() 是ggplot2中一种作图的软件。作用是通过映射,用条形图的高度来展示数值多少
    #aes(x=category,fill=strand)    表示映射关系,横轴展示分类,通过填充不同颜色来区分strand 这里也可以吧映射放到ggplot(fp,aes())里面,但是这样相当于全局映射,不便于添加其他元素和映射
    #position="fill"    表示strand里面两个元素的展示方式。可以有三种方式"stack"(堆叠),"dodge"(并列显示)和"fill"(百分比填充),默认为"stack"。
    #这里面是ggplot自动统计来正负链个数
    

    5.5 R语言绘制富集分析柱形图

    相关文章

      网友评论

        本文标题:学习生物信息,我都应该掌握那些内容?

        本文链接:https://www.haomeiwen.com/subject/kzblpqtx.html