一. Myrna的简介
Myrna主要用途是,使用RNA-seq数据来计算基因差异表达。
Myrna的运行步骤是:
1. 将reads进行预处理; 2. 使用Bowtie将reads比对到参考序列上; 3. 使用R/Bioconductor将比对结果定位到genes上; 4. 标准化基因的比对数; 5. 计算差异表达的统计数据; 6. 将结果绘制出图。
Myrna有三种使用方式:
1. 在云端使用,Amazon's Elastic MapReduce service; 2. 在 hadoop 集群上使用; 3. 在单独的一台计算机上使用。
参考文献(2010):Cloud-scale RNA-sequencing differential expression analysis with Myrna
只有读该文献才能全面了解Myrna.
二. Myrna的安装
1. 下载Myrna并安装
cd /home/chenlianfu/programs/ $ wget http://nchc.dl.sourceforge.net/project/bowtie-bio/myrna/1.2.0/myrna-1.2.0.zip $ unzip myrna-1.2.0.zip 设置Myrna的环境变量MYRNA_HOME: $ vim ~/.bashrc export MYRNA_HOME=/home/chenlianfu/programs/myrna-1.2.0
2. 安装bowtie
安装完bowtie,再将其执行程序bin目录添加到环境变量MYRNA_BOWTIE_HOME。或者在运行程序时候使用参数 –bowtie
设置Myrna的环境变量MYRNA_BOWTIE_HOME: $ vim ~/.bashrc export MYRNA_BOWTIE_HOME=/home/chenlianfu/programs/bowtie-1.0.0/
3. 安装R,Bioconductor及必须的包
本软件(版本1.20)的使用必须使用R2.14,高版本的R还不兼容。故使用Myrna自带的一个脚本来下载R2.14并安装。
# R 安装R包 lmtest 和 multicore: > install.packages("lmtest") > install.packages("multicore") 安装Bioconductor: > source("http://bioconductor.org/biocLite.R") > biocLite() 安装Bioconductor包 IRanges,geneplotter,BiocGenerics 和 biomaRt: > source("http://bioconductor.org/biocLite.R") > biocLite("IRanges") > biocLite("geneplotter") > biocLite("BiocGenerics") > biocLite("biomaRt") 设置Myrna的环境变量MYRNA_RHOME(包含bin/R和bin/Rscript): vim ~/.bashrc export MYRNA_RHOME=/usr/local/
4. 安装SRA toolkit
如果输入Myrna的软件是.sra格式,则需要SRA toolkit中的程序 fastq-dump 了。安装玩 SRA toolkit 后,再设置Myrna的环境变量:
$ vim ~/.bashrc export MYRNA_SRATOOLKIT_HOME=/home/chenlianfu/programs/sratoolkit.2.3.1-centos_linux64/bin/
5. 测试本地化运行程序
使用本地化运行程序myrna_local来进行测试,通过表示安装成功。
$ $MYRNA_HOME/myrna_local --test ... PASSED install test
6. 设置的一些环境变量
$ vim ~/.bashrc export MYRNA_BOWTIE_HOME=/home/chenlianfu/programs/bowtie-1.0.0/ export MYRNA_HOME=/home/chenlianfu/programs/myrna-1.2.0 export MYRNA_RHOME=/usr/local/ export MYRNA_SRATOOLKIT_HOME=/home/chenlianfu/programs/sratoolkit.2.3.1-centos_linux64/bin/ export MYRNA_REFS=/home/chenlianfu/programs/myrna-1.2.0/reftools/
三. Myrna的使用
1. 本地化运行程序的使用参数:
--reference <path> 参考序列的文件夹。该文件夹由一个参考序列jar文件使用unzip解压缩出来的,里面主 要包括两个子文件夹:1."index":里面包含bowtie-bulid产生的index文件;2."ival s":里面包含gene的注释信息。 --input <path> 输入的文件或文件夹路径。如果不指定--preprocess 参数,则该参数后接的是一个 文件夹,该文件夹包含了reads的预处理结果,即 --preprocess-output 的输出结果; 如果指定了 --preprocess 参数,则是一个文本文件,内容为输入的序列文件的一个清单 (Labeled manifest files)。该文件描述了输入序列文件和样本的信息。序列文件是 FASTQ或sra格式文件,FASTQ格式可以是gzip或bzip2压缩。该文件每一行代表一个测序 的结果,是以下是 unpaired reads 和 paired reads 的书写方法: <URL>(tab)<Optional MD5>(tab)<Sample label> <URL 1>(tab)<Optional MD5 1>(tab)<URL 2>(tab)<Optional MD5 2>(tab)<Sample label> 其中URL是测序结果文件的路径,可以是ftp站点上的文件;MD5是可选的,以利于下载文件后 用于文件的完整度检验,若不想进行检验,则设置为0;<Sample label>的写法为: <Group ID>-<BioRep ID>-<TechRep ID> 其中,Group ID 是样品的名称;BioRep ID是生物学重复的代号;TechRep ID是技术性 重复的代号。比如: Tumor-Subject1-Lane1 Tumor-Subject1-Lane2 Tumor-Subject2-Lane1 Tumor-Subject2-Lane2 Normal-Subject1-Lane1 Normal-Subject1-Lane2 Normal-Subject2-Lane1 Normal-Subject2-Lane2 --output <path> 输出文件的路径。如果 --just-preprocess 参数给定,则只是对reads文件进行预 处理,给出结果。默认情况下,则是输出最终结果:一系列的表格,图和比对文件。 --intermediate <path> default: /tmp/myrna/intermediate/<PID> 中间结果文件存放的临时路径。如果指定了 --keep-intermediates 或 --keep- all参数,则该文件则会一直存在。默认情况下会使用完毕后删除。 --preprocess-output <path> 对reads进行预处理的输出路径。由于对reads的预处理会消耗很多时间,因此保留数据 的预处理结果,对再次运行程序来说,可以节约很多时间。 --keep-intermediates <path> 保留程序运行中生成的结果文件夹和文件,默认情况下这些文件会被尽快的自动删除。 --keep-all 保留所有的临时文件,默认情况下这些文件会被尽快的自动删除。 --cpus <int> default: 1 设定程序运行使用的CPU线程数。 --bowtie <path> 设定Bowtie的路径,以运用bowtie进行比对的步骤,该参数优先于环境变量MYRNA_ BOWTIE_HOME,$MYRNA_HOME/bin的子目录,以及 PATH。 --fastq-dump <path> 设定 SRA toolkit 中的 fastq-dump 程序的位置,优先于PAHT。 --Rhome <path> 设定R和R/Bioconductor安装程序的位置,该目录包括 bin/R 和 bin/Rscript, 优先于环境变量MYRNA_RHOME和PATH。
2. emr, hadhoop 和 local运行的共同参数
--quality { phred33 | phred64 | solexa64 } default: phred33 设置输入reads的碱基质量格式。phred33:输入的碱基质量等于ASCII码值加上33;最 低碱基质量是“#”. phred64:输入的碱基质量等于ASCII码值加上64;最低碱基质量是“B”. --preprocess 当 --input 的参数是一个文本文件(含有输入序列信息的清单)时,则必须加入该参数 ;该参数指定出了文本文件中的序列信息来,Myrna则进行这些reads的预处理,将生成的预处 理reads结果输出到 --preprocess-output 指定的文件夹;如果不指定--preprocess -output参数,则输出到 --intermediate 指定的文件夹。 --just-preprocess 和 --preprocess 参数类似,但不同的是进行reads预处理后,即停止pipeline的 运行,并将结果输出到 --output指定的文件夹。 --just-align 不将Myrna的pipeline运行到底,只是运行到比对结束,并将结果存放到 --output 参数指定的URL. --resume-align 使用此命令来运行上一个命令(比对完毕)过后的阶段。 此时 --input 的参数则是上 一个命令中 --output的URL. --bowtie-args "<args>" default: -m 1 设置bowtie的参数。 --discard-reads <fraction> default: 0.0 丢弃一定比例输的reads来进行程序的运行。比如 0.5 则代表丢弃50%的reads。该参 数应用于所有输入的reads。对于程序的调试比较有用。 --influence <int> default: 1 一个read的"influence"长度。和下面3个参数有关系。 --from3prime 设定修剪过后的reads从3'端开始衡量其"influence"。此参数是默认设置,和下2个 参数是冲突的,三选一。 --from5prime 设定修剪过后的reads从5'端开始衡量其"influence"。 --from-middle 设定修剪过后的reads从中间开始衡量其"influence"。 --top <int> default:50 将genes按p值从小到大排序,Myrna将报告出前一定数目的genes的比对结果。默认是 50个genes。 --family { poisson | gaussian } default: poisson 设定检测gene的差异表达的统计方法。possion适合于样本<=10个;gaussian适合于 样本数>10。当然这不是硬性规定,具体还要用实验验证。当 --nulls 参数是 0 或者 没有 给出时,则使用此参数中的模型来作参数检验(parameric statistical model),计算 出p-value.当 --nulls参数大于0,则使用permutation test来计算p-value. --normalize { total | median | upper-quartile } Each count (or pooled count) is "normalized" according to a per-label baseline normalization factor. --gene-footprint { union | intersect } default: intersect 计算比对和genges的关联时所使用的 gene footprint。intersect 和 union 分别对应这reference文件夹中的ivals子文件夹中的ui和un文件夹;其中un对应的基因区 要比ui多,并且un的基因区包含了ui的,个人推测可能是un包含了5'端和3'端非翻译区。 union: 将一个比对结果指向某一个gene的必要条件是,read的"influence"(默认下是 从3'端开始)有一个碱基和一个gene任意的(any)transcripts重叠;intersect: 将 一个比对结果指向某一个gene的必要条件是,read的"influence"(默认下是从3'端)有一 个碱基和一个gene所有的(all)的transcripts重叠。我个人的理解,通俗来讲:以上数个 参数的目的是计算比对到gene上的reads的数目。根据比对结果,默认将reads从3'端开始 计算,如果比对上reads有一个碱基和一个gene的一个外显子重叠,则在union模式下这个 gene则可以得到一个read比对结果;而在intersect模式下,这个reads必须和该gene所 有的外显子有一个碱基的重叠,则该gene才能得到一个read比对结果。可见intersect要求 read贯穿整个gene才行能得到比对到gene的数目;因此,对于short reads进行基因表达 量计算的时候,则需要选union模式,不能选默认的模式。 --paired-ttest 如果设定此值,则表示: 1, 是对2个样作比较; 2, 这2个样的生物学重复必须要匹 配。当该参数使用时,将对两个样作t检验。 --nulls <int> Set the number of null statistics generated per gene to <int>. A null statistic for a gene is calculated by permuting the sample labels for the counts and normalization factors, then re-calcula- -ting the statistical test. 如果<int>是一个非0的数,则p-values由permu- -tation test得出;否则,p-values由参数检验模型(parametric statistical model)计算出。 --truncate <int> 将所有的reads都从3'端开始修剪,直到该参数设定的长度;小于该长度的reads依然 保留。 --truncate-discard <int> 将所有的reads都从3'端开始修剪,直到该参数设定的长度;而小于该长度的reads全 部剔除。 --ditch-alignments 不输出比对结果和图形结果;依然报告差异表达基因,p-values和相关的图形结果。 --discard-mate <int> 对于 paired-end reads,去掉其中一端的序列( 1 或 2 )再来进行比对。这对于 一些数据集的比较有用,比如数据集中同时有 paired reads 和 unpaired reads时。 --pool-reps 将所有的重复,包括生物学重复和技术性重复融合后,再进行计算统计。 --pool-tech-reps 将技术性重复进行融合后,在进行计算统计。 --test 搜索Myrna需要的工具(Bowtie, R/Bioconductor 和 fastq-dump),来检测 Myrna能否正确运行。 --tmpdir `<int>` default: /tmp/Myrna/invoke.scripts 临时文件存放的位置,是一个动态生成的脚本。
四. 个人感受
使用时,生成的结果没有rpkm值,p-value值及预期的图表结果。不知到哪儿出错了。比对的之前的过程都是正常的。同时使用sra文件时候,提示fastq-dump找不到或无法执行(实际上在–test的时候可行)。无法正常使用myrna。费心费力,不想深究了…悲剧…再去研究cufflinks去。