一:下载安装该软件
软件:NCBI提供的igblastn(linux环境)
需要自己去NCBI的ftp里面下载
ftp://ftp.ncbi.nlm.nih.gov/blast/executables/igblast/release/
要保证igblastn程序文件和以下三个文件夹在同一目录,可以自行下载ncbi的igblast程序,同时要下载这些东西。 Continue reading
一:下载安装该软件
软件:NCBI提供的igblastn(linux环境)
需要自己去NCBI的ftp里面下载
ftp://ftp.ncbi.nlm.nih.gov/blast/executables/igblast/release/
要保证igblastn程序文件和以下三个文件夹在同一目录,可以自行下载ncbi的igblast程序,同时要下载这些东西。 Continue reading
一:下载安装该软件
下载安装该软件: wget https://codeload.github.com/TransDecoder/TransDecoder/tar.gz/2.0.1
解压进入该目录,查看里面的文件
make一下就可以用了,看起来好像是依赖于perl模块的
这个TransDecoder.LongOrfs就是我们这次需要的程序,查看该程序,的确真是一个perl程序,看来perl还是蛮有用的。
二:准备数据
它里面有个测试数据,是比较全面的,也比较复杂,我就不贴出来了,反正我是那trinity组装好的fasta格式的转录组数据来预测ORF的。
三:运行命令
它给的测试命令也很复杂
## generate alignment gff3 formatted output
../util/cufflinks_gtf_to_alignment_gff3.pl transcripts.gtf > transcripts.gff3
## generate transcripts fasta file
../util/cufflinks_gtf_genome_to_cdna_fasta.pl transcripts.gtf test.genome.fasta > transcripts.fasta
## Extract the long ORFs
../TransDecoder.LongOrfs -t transcripts.fasta
当然我们只需要看最后一步,这是重点
我这里是直接对我们的trinity组装好的转录本进行预测ORF
/home/jmzeng/bio-soft/TransDecoder/TransDecoder.LongOrfs -t Trinity.fasta
命令很简单
输出来的文件就有预测的蛋白文件,这个文件是trinotate对转录本进行注释所必须的文件
四:输出文件解读
longest_orfs.cds 这个是预测到的cds碱基序列,
longest_orfs.gff3 这个是预测得到的gff文件
longest_orfs.pep 这个就是预测得到的蛋白文件
一,下载该软件
wget http://ftp-trace.ncbi.nlm.nih.gov/sra/sdk/current/sratoolkit.current-ubuntu64.tar.gz
tar xzf sratoolkit.current-centos_linux64.tar.gz
解压直接使用即可,里面有一大堆的软件,针对不同的测序仪,不同的数据 Continue reading
查某个基因家族在某物种的具体信息
我很伤心,不知道是不是我写的教程还是不够人性化,一个朋友在群里面问如何知道NAC基因家族在拟南芥里面的105个基因信息,我随便给他示范了一下在人类里面如何找,希望他能触类旁通,结果他不会linux,啥生信基础都没有,我只会诱导他简单学习一下,希望他至少明白什么的taxid。所以我给了他我之前写的教程,只希望他告诉我拟南芥的taxid我就帮他把那105个基因找出来。 Continue reading
生信常用论坛seq-answer里面所有帖子爬取
这个是爬虫专题第二集,主要讲如何分析seq-answer这个网站并爬去所有的帖子列表,及标签列表等等,前提是读者必须掌握perl,然后学习perl的LWP模块,可以考虑打印那本书读读,挺有用的!
其实爬虫是个人兴趣啦,跟这个网站没多少关系,本来一个个下载,傻瓜式的重复也能达到目的。我只是觉得这样很有技术范,哈哈,如何大家不想做傻瓜式的操作可以自己学习学习,如果不懂也可以问问我!
http://seqanswers.com/这个是主页
http://seqanswers.com/forums/forumdisplay.php?f=18 这个共570个页面需要爬取
其中f=18 代表我们要爬去的bioinformatics板块里面的内容
http://seqanswers.com/forums/forumdisplay.php?f=18&order=desc&page=1
http://seqanswers.com/forums/forumdisplay.php?f=18&order=desc&page=570
<tbody id="threadbits_forum_18">这个里面包围这很多<tr>对,
前五个<tr>对可以跳过,里面的内容不需要
生信常用论坛bio-star里面所有帖子爬取
这个是爬虫专题第一集,主要讲如何分析bio-star这个网站并爬去所有的帖子列表,及标签列表等等,前提是读者必须掌握perl,然后学习perl的LWP模块,可以考虑打印那本书读读,挺有用的!
http://seqanswers.com/ 这个是首页
http://seqanswers.com/forums/forumdisplay.php?f=18 这个共570个页面需要爬取
http://seqanswers.com/forums/forumdisplay.php?f=18&order=desc&page=1
http://seqanswers.com/forums/forumdisplay.php?f=18&order=desc&page=570
<tbody id="threadbits_forum_18">这个里面包围这很多<tr>对,
前五个<tr>对可以跳过,里面的内容不需要
找橡胶测序数据无果
所以我只好找了他们所参考的草莓(strawberry, Fragaria vesca (2n = 2x = 14),a small genome (240 Mb),)的文章,是发表是nature genetics上面的
http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3326587/
可以看到它的SRA索取号。
草莓组装结果:Over 3,200 scaffolds were assembled with an N50 of 1.3 Mb .
Over 95% (209.8 Mb) of the total sequence is represented in 272 scaffolds.
草莓基因息:Gene prediction modeling identified 34,809 genes, with most being supported by transcriptome mapping.
草莓染色体信息:Paradoxically, the small basic (x = 7) genome size of the strawberry genus, ~240 Mb,
offers substantial advantages for genomic research.
草莓来源:diploid strawberry F. vesca ssp. vesca accession Hawaii 4
(National Clonal Germplasm Repository accession # PI551572).
然后我去NCBI上面下载这三个数据
SRA020125 共有四个数据:
http://www.ncbi.nlm.nih.gov/sra/SRX030575[accn] | Total: 4 runs, 4.7M spots, 2.6G bases, 5.5Gb |
http://www.ncbi.nlm.nih.gov/sra/SRX030576[accn] (3 KB PE) | Total: 2 runs, 2.2M spots, 908.5M bases, 2.1Gb |
http://www.ncbi.nlm.nih.gov/sra/SRX030577[accn] (20KB片段) | Total: 2 runs, 1.9M spots, 800M bases, 1.8Gb |
http://www.ncbi.nlm.nih.gov/sra/SRX030578[accn] | Total: 3 runs, 4M spots, 2.2G bases, 4.6Gb |
挂在后台自动下载
好了,有了这些数据我们就要进行基因组的一系列分析啦!!!
不过我们可以先看看他们这个研究小组的成果
首先他们建造了一个关于草莓的基因组信息网站
https://strawberry.plantandfood.co.nz/
跟我之前在水科院做鲫鱼鲤鱼的差不多
直接在里面就可以下载他们做好的所有数据,也可以可视化。
它的染色体如下,非常简单,就七条染色体
http://www.rosaceae.org/species/fragaria/fragaria_vesca/genome_v1.1
我找到了它组装好的草莓基因组地址,用批处理全部下载了
研读橡胶的基因组文章
我本科的前两年在海南儋州读书,那时候旁边就是橡胶所,很多同学也在那边做毕业论文什么的,我一直以为那里是全世界的橡胶中心,所有的先进技术都在那里产生,结果,前些天跟一个橡胶所的老师聊天才发现,居然橡胶(Hevea brasiliensis)的基因组已经发表了,可是,跟橡胶所没有半毛钱关系,更搞笑的事情是,堂堂一个基因组文章居然发表在BMC这样的杂志,真不知道是基因组的年代已经过去了还是他们做的实在是太差了,反正我看不过去了,所以研读他们的文章,并且下载数据测试一下。
文章地址如下:http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3575267/
可以看到它过于数据的描述都在补充材料1里面,所以我下载了补充材料。
可以看到所有的测序数据的描述,45个G的i llumina的200bp的双端测序,27个G的illumina的200bp的双端测序,约10G左右的长片段(8kb,20kb)罗氏454数据,最后还有一点点solid数据,它这样的测序策略好像是模仿的2011年发布的草莓基因组数据。
但是补充材料里面没有列出下载地址,我有点困惑!
按照道理我研读文献的步骤应该没有错,有可能是因为这个文章发表的杂志水平太低,所以不要求他们把测序原始数据上传到NCBI的SRA里面。或者是他们本身觉得文章发的不够档次,不想公布数据,所以先留着自己做精细分析,等发了大文章再公布原始数据。
然后我在NCBI的SRA里面查找了关于橡胶的原始数据,果真没有
仅有的10个数据,都是别的小组做的RNA-seq的内容。
De novo transcriptome analysis of abiotic stress responsive transcripts of Hevea brasiliensis.
所以我只好找了他们所参考的草莓(strawberry, Fragaria vesca (2n = 2x = 14),a small genome (240 Mb),)的文章,是发表是nature genetics上面的
我真的不是程序员,也没时间去自己写一个APP,无意中看到了一个APP的弹出页面写着简易APP工厂支持,我试着搜索了一下,才知道,原来他们提供了一个平台,傻瓜式的创建一个自己的APP,当然,现在好像只是免费提供安卓版本,不过也非常实用!!!
非常easy的节目,大家如果有兴趣,也可以自己下载一个!!!
然后我顺便搜索了一下我的网站效果,发现现在终于被百度搜录了,而且,居然,我很久以前写的菜鸟生物信息学居然还能排名第二,我很久以前的想法就是分享一下自己学习过程中的艰辛曲折,给后学者们借鉴,希望这样可以帮到更多的朋友!
希望有在深圳的生信从业人员或者学生能看到此广播,我们可以组成兴趣小组交流一下各自所学,或者合作翻译一些技术文档或者制作生信常用软件的使用说明书。
YGC是我非常喜欢的一个博客,但是之前都是英文,所以我没怎么看,余老师还在R领域很有建树,个人发表了5个R在生信方面的包,我后面有空会一个个试用学习。
本文转自http://ygc.name/2014/08/27/insertion-size/
里面详细解释了NGS测序的几个基本概念,也是我之前一直弄混淆的概念,包括插入片段、单端测序,双端测序,配对测序,contig,scaffold等等
在进行测序的时候,需要将DNA打断,构建library,这些fragment需要接上adaptor,好进行扩增,illumina的测序,可以有single end和paired end两种,分别从一端和两端进行测序。
fragment ======================================== fragment + adaptors ~~~========================================~~~ SE read ---------> PE reads R1---------> <---------R2 unknown gap ....................
insertion并不是指R1和R2之间的unknown gap,早在NGS之前,当我们在使用ecoli构建载体的时候,这个概念就已经形成,它是adaptors之间的序列。而unknown gap则称之为inner mate:
PE reads R1---------> <---------R2 fragment ~~~========================================~~~ insert ======================================== inner mate ....................
显然我们不希望看到大量的unknown gap,所以要制造短的fragment,而且技术不断发展,测序长度也越来越长,于是可以测通fragment:
fragment ~~~========================================~~~ insert ======================================== R1 -------------------------> R2 <----------------------- overlap :::::::::: stitched SE read --------------------------------------->
这样R1和R2就有overlap,合并一致序列,就可以得到完整的fragment,使用短的fragment,也就是insertion size比较小的library,测序的结果coverage比较大,因为我们可以测通fragment.
虽然adaptor不会被测序,但如果fragment太短,被读通了,则另一端的adaptor就会被测到。
tiny fragment ~~~~========================~~~~ insert ======================== R1 --------------------------> R2 <-------------------------- read-through !!! !!!
如果MiSeq设置正确的话,读通的adaptor是会被切除了,这样就会获得长度不一致的short reads,也可以使用N来替换adaptor序列,这样长度一样,但会在5' end看到很多N。如果没设置好,reads里含有adaptor序列,那么必须要通过软件去除,否则后续的分析都会有问题。
所以insertion size小有个好处,测序的genome coverage高,但是在进行de novo assembly的时候,有一个问题,如果基因组含有比read length还要长的重复元件时,就无法拼接,所以得到的是很多的contigs,它们之间的gap要长于insertion size且无法确定。这个问题是相当普遍的,即使是相对简单的ecoli基因组,也有一定数量的重复元件。
这个问题需要使用大的insertion size进行paired end测序来解决。
fragment + adaptors ~~~========================================~~~ PE reads R1---------> <---------R2 unknown gap ....................
在这种insertion size比较大的情况下,我们可以估计R1和R2之间的距离,只要有一个片段能够被mapped到unique position的话,那么另一个片段的大致位置就可以确定。所以为了达到好的拼接效果,长fragment的library也是必须的,它有可能给出 contigs间的相对位置。
所以理想的情况是使用multiple insert libraries,short-insert library可以保证获得足够的coverage,它可以告诉你contigs之间的序列,但信息是local的,它没办法告诉你怎么拼;而long- insert libary则可以告诉你一些相对global的信息。
在上面这个测试的数据中,加了long-insert libary虽然在coverage上没多少变化,但N50和最大的contig都显著提高,4.5Mb已经覆盖了~98%的ecoli基因组。
- See more at: http://ygc.name/2014/08/27/insertion-size/#sthash.fKmWKoaf.dpuf
转录组cummeRbund操作笔记
这是跟tophat和cufflinks套装紧密搭配使用的一个R包,能出大部分文章要求的标准化图片。
一:安装并加装该R包
安装就用source("http://bioconductor.org/biocLite.R") ;biocLite("cummeRbund")即可,如果安装失败,就需要自己下载源码包,然后安装R模块。
然后把cuffdiff输出的文件目录拷贝到R的工作目录,或者自己设置工作目录
二:读取FN目录下面的所有文件。
可以看到把cuffdiff下面的文件夹所有的文件都读取到了,里面有如下文件,包括genes,isoforms,cds,tss这四种差异情况都读取了。
三:表达水平分布图
csBoxplot(genes(cuff_data))
画出热图如下
六、得到差异的genes,isoforms,TSS,CDS等等
diffData <- diffData(myGenes )
只有一百个有表达差异的基因
最后贴出一个综合性的代码,算了,太浪费空间了,把整个空间搞得不好看,就不贴了。
这个代码可以自动运行出图;
转录组edgeR分析差异基因
edgeR是一个研究重复计数数据差异表达的Bioconductor软件包。一个过度离散的泊松模型被用于说明生物学可变性和技术可变性。经验贝叶斯方法被用于减轻跨转录本的过度离散程度,改进了推断的可靠性。该方法甚至能够用最小重复水平使用,只要至少一个表型或实验条件是重复的。该软件可能具有测序数据之外的其他应用,例如蛋白质组多肽计数数据。可用性:程序包在遵循LGPL许可证下可以从Bioconductor网站。
一:下载安装该软件
下载安装edgeR这个R包,因为这是一次讲R包的下载,我就啰嗦一点,这种生物信息学的包不同于普通的R包,是需要用biocLite来安装的,命令如下
转录组HTseq对基因表达量进行计数
一:下载安装该软件
下载htseq这个python模块安装解压包,依赖于很多python的其它安装包及库,模块,我最讨厌python了,在有些电脑上特别难安装,而且服务器还有权限的问题。
解压进入该目录,输入 python setup.py install --user 记住,是- - 而不是—
这样只是把这个软件安装到自己的目录
安装完毕后,会出现这两个程序,在自己的python库里面,可以直接调用这两个程序的,我这里它们的路径是 .local/bin ,很奇怪的一个路径,我也是用find命令才找到的
仿写fastqc软件的一些功能(下)
文件来自于上面perl代码的输出文件,好像算法有点问题,26G的文件居然处理近一个小时才出数据!
R语言本身自带的画图工具都很丑,懒得说了,可以用ggplot2来重新画一个,不是项目要求没有报酬我就懒得画了,大家面前看看画图原理即可。
转录组cufflinks套装的使用
cufflinks套装有很多,我们主要使用的只有三个
Cufflinks是用来处理tophat的输出的bam文件然后输出gtf文件
cuffmerge把多个样本的gtf文件合并的,也没啥子用,主要是测多个样本可能会需要
cuffdiff算出分组的bam文件里面的差异基因。
一:下载安装该软件
是二进制版本,找到网址,然后用wget下载,解压即可使用
转录组比对软件tophat的使用
为什么要用这个软件?:因为转录组reads比对到基因组reads用bwa和bowtie的效果都不够好,所以我们选择tophat
它做了什么?:tophat把测序的转录组的原始reads比对到了参考基因组上面,并且输出了bam(二进制的sam)文件比对结果给我们。(fastq--->bam)
一:下载安装该软件
其实一般的生信服务器自然会有高手给安装好了,你只需调用即可,这里我给大家演示一下如何安装。
wget http://ccb.jhu.edu/software/tophat/downloads/tophat-2.0.13.Linux_x86_64.tar.gz
仿写fastqc软件的部分功能(上)
前面我们介绍了fastqc这个软件的使用方法 http://www.bio-info-trainee.com/?p=95 ,这是一个java软件,但是有些人服务器没有配置好这个java环境,导致无法使用,这里我贴出几个perl代码,也能实现fastqc的部分功能
统一测试文件是illumina的phred33格式的fastq文件,共100000/4=25000条reads,读长都是101个碱基
程序名-fastq2quality.pl
使用命令:perl fastq2quality.pl SRR504517_1.fastq >quality.txt
功能: 把fastq格式的每条原始reads的第四行ascii码质量值,转换为Q值并输出一个矩阵,有多少条reads就有多少行,每条reads的碱基数就是列数。
用R语言画3D条形图
我暂时找到的能画3D条形图的R语言包就有四个,接下来我们分别比较一下它们的优劣
library(rgl)
library(scatterplot3d)
library(epade)
library(latticeExtra)
数据塑形包(reshape2)也需要加载
生信菜鸟养成手册
背景:生物小本,懂做一些分子实验,了解一些生物背景知识。
目标:成为生信菜鸟,找到一份生物信息学相关的工作。
1、计算机基础(linux+perl+R 或者 python+matlab)
2、生信基础知识(测序+数据库+数据格式)
3、生信研究领域(全基因组,全转录组,全外显子组,捕获目标区域测序)
4、生信应用领域(肿瘤筛查,产前诊断,流行病学,个性化医疗)
Trinity进行转录组组装的使用说明
一:下载安装该软件
去官网下载trinity并解压安装 http://trinityrnaseq.github.io/
安装非常简单,一个make即可
这个软件比较大,约150M。所以安装需要一会时间,以下是安装进程日志,可以看出trinity这个软件安装的同时还附带着好几个测序一起安装进来了。