26

阅读捕获测序文章并下载数据

一.阅读文献找到SRP

该文献讲了单分子测序在医疗领域的一个应用,我感觉挺重要的,就分析了一下,然后下载了数据,准备处理一下。

Single-step capture and sequencing of natural DNA for detection of BRCA1 mutations

 

捕获测序文章解析并下载数据162

 

在NCBI查到该数据地址,并且用脚本下载即可

http://www.ncbi.nlm.nih.gov/sra/?term=SRP007097

捕获测序文章解析并下载数据348

下载之后的数据如下,共19个测序文件,都是200K左右大小,那两个一百多M的可能是下载错了

 

 

for i in {32..52}

do

wget ftp://ftp-trace.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByRun/sra/SRR/SRR258/SRR2588$i/SRR2588$i.sra

Done

 

下载的19个数据,都是只有1万多条序列。

捕获测序文章解析并下载数据346

因为这些判断都是对BRCA1这个基因进行目标性测序,所以接下来需要对它们进行特殊的处理。

26

一个基因的生信之旅

感觉大家对很多生物信息学的术语都不甚了解,我这里简单的从一个基因开始,扩展开来讲一讲生信数据库,及它相关的一些术语!

我要讲的基因是BRCA1,这是一个与乳腺癌以及卵巢癌都息息相关的基因。而BRCA1是它的英文缩写简称,也是通常学者们进行交流十它的名字。它的全称是breast cancer 1,每个基因都会有一个简称,比如下面这些,在human里面这些简称多大47732个,正常人都不会认识它们所有,只需要碰到了去数据库搜索即可,但是搞医疗健康的,必须熟悉癌症50基因。

一个基因的生信之旅247

这样的缩写简称其实弊端很多,单词毕竟是有限的,而且缩写也没有语义。所以NCBI给每个基因都定义了一个entrez ID号,是整数的排序,具体大家可以去看NCBI发的一篇文献,专门讲解了entrez ID号的好处。

1 A1BG

2 A2M

3 A2MP1

9 NAT1

10 NAT2

11 NATP

12 SERPINA3

13 AADAC

14 AAMP

这里我们来找一下我们的BRCA1这个基因在生物信息数据库里面的其它信息,在NCBI的ftp里面有一个文件是Homo_sapiens.gene_info里面包含着人类所以基因的全部信息

9606  首先这个基因在human上面的,而human被NCBI定义的taxid是9606

672  然后这个基因的被NCBI定义的entrez ID号是672

BRCA1  这个当然就是这个基因的英文缩写名称啦

-      这个表明这个基因在负链什么

BRCAI|BRCC1|BROVCA1|FANCS|IRIS|PNCA4|PPP1R53|PSCP|RNF53

这个可能是基因以前的名称,或者是在其它研究领域的一些名称。MIM:113705|HGNC:HGNC:1100|Ensembl:ENSG00000012048|HPRD:00218|Vega:OTTHUMG00000157426

这里面包含在它在其它数据库的信息,我们的NCBI用entrez ID号672来标识它,相应的ensembl数据用ensembl ID号ENSG00000012048来标识它,还有什么MIM数据库,HGNC数据库,Vega数据库我就不详细讲啦

17 17q21 这个说明它在human的17号染色体的位置信息

下面一堆都是这个基因的描述,它的功能等等。

breast cancer 1, early onset protein-coding BRCA1 breast cancer 1, early onset

O BRCA1/BRCA2-containing complex, subunit 1|Fanconi anemia, complementation group S|RING finger protein 53|breast and ovarian cancer susceptibility protein 1|breast and ovarian cancer sususceptibility protein 1|breast cancer type 1 susceptibility protein|protein phosphatase 1, regulatory subunit 53

20150201

这样我们就把好几个数据库给串起来了,也大致了解了一个基因的各种信息,但是,这样肯定是不够的。

接下来我们就不用BRCA1来称呼这个基因了,我们统一用NCBI定义entrez ID号672来称呼这个基因,当然用ensembl ID号ENSG00000012048也可以,它们都是比较通用的。

ENSG00000012048 672 这个基因在GO数据库里面可以找到67个功能信息,分别是以下

GO:0000151 GO:0000724 GO:0000724 GO:0000794 GO:0003677 GO:0003684 GO:0003713 GO:0003723 GO:0004842 GO:0005515 GO:0005634 GO:0005654 GO:0005694 GO:0005737 GO:0005886 GO:0006260 GO:0006281 GO:0006301 GO:0006302 GO:0006302 GO:0006349 GO:0006357 GO:0006359 GO:0006633 GO:0006915 GO:0006974 GO:0006978 GO:0007059 GO:0007098 GO:0008270 GO:0008274 GO:0008630 GO:0009048 GO:0010212 GO:0010575 GO:0010628 GO:0015631 GO:0016567 GO:0016874 GO:0019899 GO:0030521 GO:0030529 GO:0031398 GO:0031436 GO:0031572 GO:0031625 GO:0035066 GO:0035067 GO:0042127 GO:0042981 GO:0043009 GO:0043234 GO:0043627 GO:0044030 GO:0044212 GO:0045717 GO:0045739 GO:0045766 GO:0045892 GO:0045893 GO:0045893 GO:0045944 GO:0045944 GO:0046600 GO:0050681 GO:0051571 GO:0051572 GO:0051573 GO:0051574 GO:0051865 GO:0070512 GO:0070531 GO:0071158 GO:0071356 GO:0071681 GO:0085020 GO:1902042 GO:2000378 GO:2000617 GO:2000620

由于GO太多了,我简单讲几个

ubiquitin ligase complex

double-strand break repair via homologous recombination

double-strand break repair via homologous recombination

condensed nuclear chromosome

DNA binding

damaged DNA binding

transcription coactivator activity

RNA binding

ubiquitin-protein transferase activity

protein binding

都是描述这个基因的功能的。

到这里我们大致了解了这个基因的功能,但是还不够。

然后可以查到它有一下6个转录本,都有二十多个外显子。

NR_027676

NM_007300

NM_007299

NM_007298

NM_007297

NM_007294

在hg19这个参考基因组的起始终止坐标,还有各个外显子的起始终止坐标都能找到。

41196311,41199659,41201137,41203079,41209068,41215349,41215890,41219624,41222944,41226347,41228504,41234420,41242960,41243451,41247862,41249260,41251791,41256138,41256884,41258494,41267742,41276033,41277198

 

41197819,41199720,41201211,41203134,41209152,41215390,41215968,41219712,41223255,41226538,41228631,41234592,41243049,41246877,41247939,41249306,41251894,41256278,41256973,41258550,41267796,41276132,41277340

http://asia.ensembl.org/Homo_sapiens/Gene/Summary?db=core;g=ENSG00000012048;r=17:43044295-43125483

在ensembl里面关于这个基因的描述如下。

breast cancer 1, early onset [Source:HGNC Symbol;Acc:HGNC:1100]

BRCC1, FANCS, PPP1R53, RNF53

Chromosome 17: 43,044,295-43,125,483 reverse strand.

chromosome:GRCh38:CM000679.2:43044295:43125483:1

This gene has 29 transcripts (splice variants), 63 orthologues, is a member of 4 Ensembl protein families and is associated with 11 phenotypes.

RefSeq Gene ID 672

Uniprot identifiers: P38398

而且ensembl里面可以可视化这个基因的所有信息。

然后简单检索一下关于这个BRCA1基因的文献发表状况,居然多达2111篇文献,看来这个基因很火呀!!!

awk '{if ($1==9606 && $2==672) print }' gene2pubmed |wc

9606 672 1676470

9606 672 2001833

9606 672 2270482

9606 672 4506230

9606 672 7481765

9606 672 7545954

9606 672 7550349

9606 672 7795652

9606 672 7894491

9606 672 7894492

第三列1676470等编号是pubmed数据库的文献编号,可以直接找到关于这个基因的文献发表情况。

而直接在NCBI的pubmed数据库里面可以搜到多达11339篇文献。

esearch -db pubmed -query 'BRCA1'

Esearch这个程序是NCBI提供的,挺好用的,希望大家可以熟悉一下。

esearch -db pubmed -query 'BRCA1' | efetch -format docsum |   xtract -pattern DocumentSummary -present Author -and Title     -element Id -first "Author/Name" -element Title  >BRCA1.pubmed

一个基因的生信之旅4634

用这个代码,可以找到所有关于这个BRCA1基因的文献的作者及标题,这样可以统计在这个基因领域的研究者最出名的是谁。

至于这个基因的序列,及其转录本翻译的蛋白我就不列了,太长了,而且占位子

 

 

 

 

24

旧版本blast详解

其实我现在一般都用的是blast++了,也专门写了篇日志介绍它!

但是看到一些就的服务器上面只有blast,所以就搜了一些它的用法。

主要参考 http://www.bio.ku.dk/nuf/resources/BLAST_index.htm

很简单的两个步骤

首先建库formatdb -i Cad16_aa.fasta -p T -o F

就是把 Cad16_aa.fasta这个序列文件变成blast专用的库,-p选项中的T是代表蛋白库

然后就比对咯,比对程序有六个,需要用-p来选择

blastall -p blastx -d nr -i 19A.fa -o 19A.outm -v 1 -b 1 -m 8

上面这个命令就是选择了blastx这个比对程序,数据库是nr ,输入的查询序列是 19A.fa

然后我们输出格式的m8,这个格式很重要,我们还可以设置-a控制cpu数量,和-e控制阈值

BLAST programs

blastp Protein query > Protein database
blastn Nucleotide query > Nucleotide database
blastx Nucleotide query > Protein database (via translated query)
tblastn Protein query > Nucleotide database (via translated database)
tblastx Nucleotide query > Nucleotide database (via translated query and database) 

 

Formatting database for local BLAST

- Show a list of all arguments.
-i Input file(s) for formatting. Optional.
-p Type of file [T/F]. T = protein, F = nucleotide. Default = T.
-o Parse option [T/F]. T = Parse SeqId and create indexes, F = Do not parse or create indexes.
23

与hg19的突变相关的一些数据解释。

http://statgenpro.psychiatry.hku.hk/limx/kggseq/download/resources/

这个网站收集了大部分资料,我们就用它的,如果它倒闭了,大家再想办法去搜索吧。

与hg19的突变相关的一些数据解释210

其实这些文件都是基于NCBI以及UCSC和ensembl数据库的文件用一些脚本转换而来的,都是非常简单的脚本。

首先我们看看humandb/hg19_refGene.txt 这个文件,总共2.5万多个基因的共5万多个转录本。

与hg19的突变相关的一些数据解释325

19     可能是entrez ID,但是又不像。

NM_001291929    参考基因名

chr11   染色体

-

89057521

89223909

89059923

89223852

17      89057521,89069012,89070614,89073230,89075241,89088129,89106599,89133184,89133382,89135493,89155069,89165951,89173855,89177302,89182607,89184952,89223774,       89060044,89069113,89070683,89073339,89075361,89088211,89106660,89133247,89133547,89135710,89155150,89166024,89173883,89177400,89182692,89185063,89223909,

0

NOX4    基因的英文简称,通俗名

cmpl

cmpl

2,0,0,2,2,1,0,0,0,2,2,1,0,1,0,0,0,

然后我们看看hg19_snp141.txt这个文件

1       10229   A       -       .

1       10229   AACCCCTAACCCTAACCCTAAACCCTA     -       .

1       10231   C       A       .

1       10231   C       -       .

1       10234   C       T       .

1       10248   A       T       .

1       10250   A       C       .

1       10250   AC      -       .

1       10255   A       -       .

1       10257   A       C       .

1       10259   C       A       .

1       10291   C       T       .

1       10327   T       C       .

1       10329   ACCCCTAACCCTAACCCTAACCCT        -       .

1       10330   C       -       .

1       10390   C       -       .

1       10440   C       A       .

1       10440   C       -       .

1       10469   C       G       .

1       10492   C       T       .

1       10493   C       A       .

1       10519   G       C       .

1       10583   G       A       0.144169

1       10603   G       A       .

1       10611   C       G       0.0188246

1       10617   CGCCGTTGCAAAGGCGCGCCG   -

里面记录了以hg19为参考的所有的snp位点。

 

585

ENST00000518655 基因的ensembl ID号

chr1 + 11873 14409 14409 14409

4 基因有四个外显子

11873,12594,13402,13660, 12227,12721,13655,14409, 在基因的四个外显子的坐标

0

DDX11L1 基因的通俗英文名

none none -1,-1,-1,-1,

CTTGCCGTCAGCCTTTTCTTT·····gene的核苷酸序列

23

用annovar对snp进行注释

一、下载及安装软件

这个软件需要edu邮箱注册才能下载,可能是仅对科研高校开放吧。所以软件地址我就不列了。

用annovar对snp进行注释71

它其实是几个perl程序,比较重要的是这个人类的数据库,snp注释必须的。

用annovar对snp进行注释111

参考:http://annovar.readthedocs.org/en/latest/misc/accessory/

二,准备数据

既然是注释,那当然要有数据库啦!数据库倒是有下载地址

http://www.openbioinformatics.org/annovar/download/hg19_ALL.sites.2010_11.txt.gz

也可以用命令来下载

Perl ./annotate_variation.pl -downdb -buildver hg19 -webfrom annovar refGene humandb/

然后我们是对snp-calling流程跑出来的VCF文件进行注释,所以必须要有自己的VCF文件,VCF格式详解见本博客另一篇文章,或者搜索也行

http://vcftools.sourceforge.net/man_latest.html

三、运行的命令

首先把vcf格式文件,转换成空格分隔格式文件,自己写脚本也很好弄

perl convert2annovar.pl  -format vcf

/home/jmzeng/raw-reads/whole-exon/snp-calling/tmp1.vcf >annovar.input

用annovar对snp进行注释886

变成了空格分隔的文件

用annovar对snp进行注释900

然后把转换好的数据进行注释即可

./annotate_variation.pl -out ex1 -build hg19 example/ex1.avinput humandb/

 

四,输出文件解读

用annovar对snp进行注释1002

23

Snp-calling流程(BWA+SAMTOOLS+BCFTOOLS)

比对可以选择BWA或者bowtie,测序数据可以是单端也可以是双端,我这里简单讲一个,但是脚本都列出来了。而且我选择的是bowtie比对,然后单端数据。

首先进入hg19的目录,对它进行两个索引

samtools faidx hg19.fa

Bowtie2-build hg19.fa hg19

我这里随便从26G的测序数据里面选取了前1000行做了一个tmp.fa文件,进入tmp.fa这个文件的目录进行操作

Bowtie的使用方法详解见http://www.bio-info-trainee.com/?p=398

bowtie2 -x ../../../ref-database/hg19 -U  tmp1.fa -S tmp1.sam

samtools view -bS tmp1.sam > tmp1.bam

samtools sort tmp1.bam tmp1.sorted

samtools index tmp1.sorted.bam 

samtools mpileup -d 1000  -gSDf   ../../../ref-database/hg19.fa  tmp1.sorted.bam |bcftools view -cvNg -  >tmp1.vcf

然后就能看到我们产生的vcf变异格式文件啦!

 

当然,我们可能还需要对VCF文件进行再注释!

要看懂以上流程及命令,需要掌握BWA,bowtie,samtools,bcftools,

数据格式fasta,fastq,sam,vcf,pileup

 

如果是bwa把参考基因组索引化,然后aln得到后缀树,然后sampe对双端数据进行比对

首先bwa index 然后选择算法,进行索引。

然后aln脚本批量处理

==> bwa_aln.sh <==

while read id

do

echo $id

bwa aln hg19.fa $id >$id.sai

done <$1

然后sampe脚本批量处理

==> bwa_sampe.sh <==

while read id

do

echo $id

bwa sampe hg19.fa $id*sai $id*single >$id.sam

done <$1

然后是samtools的脚本

==> samtools.sh <==

while read id

do

echo $id

samtools view -bS $id.sam > $id.bam

samtools sort $id.bam $id.sorted

samtools index $id.sorted.bam

done <$1

然后是bcftools的脚本

==> bcftools.sh <==

while read id

do

echo $id

samtools mpileup -d 1000  -gSDf  ref.fa $id*sorted.bam |bcftools view -cvNg -  >$id.vcf

done <$1

 

 

==> mpileup.sh <==

while read id

do

echo $id

samtools mpileup -d 100000 -f hg19.fa $id*sorted.bam >$id.mpileup

done <$1

 

21

Hg19基因组的分析

下载地址我就不贴了,随便谷歌一下即可!

Genome Reference Consortium Human  ---》  GRCh3

Feb. 2009 (hg19, GRCh37)这个是重点

Mar 2006 assembly = hg18 = NCBI36.

May 2004 assembly = hg17 = NCBI35.

July 2003 assembly = hg16 = NCBI34

以前的老版本就不用看啦,现在其实都已经有hg38出来啦,GRCh38 (NCBI) and hg38(UCSC)

参考:http://age.wang.blog.163.com/blog/static/119252448201092284725460/

http://www.ncbi.nlm.nih.gov/projects/genome/assembly/grc/human/

Hg19基因组的分析570

人的hg19基因组是3G的大小,因为一个英文字符是一个字节,所以也是30亿bp的碱基。

包括22条常染色体和X,Y性染色体及M线粒体染色体。

Hg19基因组的分析643

查看该文件可以看到,里面有很多的N,这是基因组里面未知的序列,用N占位,但是觉得部分都是A.T.C.G这样的字符,大小写都有,分别代表不同的意思。

然后我用linux的命令统计了一下里面这个文件的行数,

perl -lne 'END { print $. }'  hg19.fa

awk 'END { print NR }'  hg19.fa

wc -l hg19.fa

Hg19基因组的分析834

然后我写了一个脚本统计每条染色体的长度,42秒钟完成任务!

Hg19基因组的分析1125

看来这个服务器的性能还是蛮强大的,读取文件非常快!

[perl]

while(<>){

        chomp;

        if  (/>/){

if  (exists $hash_chr{$key} ){

$len = length $hash_chr{$key};

print "$key   =>   $len\n";

}

undef %hash_chr;

$key=$_;

}

else {

$hash_chr{$key}.=$_;

}

}

[/perl]

 

然后我用seed统计了一下hg19的词频(我不知道生物信息学里面的专业描述词语是什么)

Hg19基因组的分析1171

我的程序耗费了42分钟才跑完,感觉我写的程序应该是没有问题的,让我吃惊的是总共竟然只有105万条独特的10bp短序列。然后我算了一下4的10次方,(⊙o⊙)…悲剧,原来只有1048576,之所以出现这种情况,是因为里面有N这个字符串,不仅仅是A.T.C.G四个字符。我用grep -v N seed10.txt |wc -l命令再次统计了一下,发现居然就是1048576,也就是说,任意A.T.C.G四个字符组成的10bp字符串短序列在人的基因组里面都可以找到!!!

Hg19基因组的分析1407

然后我测试了一下,还是真是这样的,真是一个蛮有意思的现象。虽然我无法解释为什么,但是根据这个结果我们可以得知连续的A或者T在人类基因组里面高频出现,而连续的G或者C却很少!

如果我们储存这个10bp字符串的同时,也储存着它们在基因组的位置,那么就可以根据这个seed来进行比对,这就是blast的原理之一!

 

 

19

转录组-GO通路富集-WEGO网站使用

一,所谓的网站,其实就是一个网页版的可视化软件接口而已

看看网站主页,看看它需要什么数据

http://wego.genomics.org.cn/cgi-bin/wego/index.pl

转录组-GO通路注释-WEGO网站使用181

二,所需要的数据

1,human.all.go.entrez,需要自己制作,每个基因名entrez ID号,对应着一堆GO通路,人有两万多个基因,所以应该有两万多行的文件。

转录组-GO通路注释-WEGO网站使用247

2,差异基因的GO通路,需要用cuffdiff得到差异基因名,然后用然后用脚本做成下面的样子。记住,上面的那个人类的背景GO文件也是一样的格式,基因名是entrez ID号,与GO通路用制表符隔开,然后每个基因所对应的GO直接用空格隔开。格式要求很准确才行。

转录组-GO通路注释-WEGO网站使用379

 

三,上传数据,出图

转录组-GO通路注释-WEGO网站使用391

点击plot画图即可,就可以出来了一个GO通路富集图

转录组-GO通路注释-WEGO网站使用420

顺便贴上wego上传数据制作的几个脚本,脚本这种东西都很难看,随便意思一下啦,用一下脚本处理就可以得到wego需要上传的数据了

1,得到差异基因名,并且转换为entrez ID号
grep yes gene_exp.diff |cut -f 3 |sort -u >diff.gene.name
cat diff.gene.name  ../Homo_sapiens.gene_info  |perl -alne '{$hash{$_}=1;print $F[1] if exists $hash{$F[2]}}' |sort -u >diff.gene.entrez
2,根据找到的差异基因的entrez ID号来找到它的GO信号,输出文件给wego网站
cat diff.gene.entrez ../gene2go |perl -alne '{$hash{$_}=1;print "$F[1]\t$F[2]" if exists $hash{$F[1]}}' |perl -alne '{$hash{$F[0]}.="$F[1] "}END{print "$_\t$hash{$_}" foreach keys %hash}' >diff.gene.entrez.go
3,得到entrez ID号跟ensembl ID号的转换hash表
perl -alne '{if (/Ensembl:(ENSG\d+)/) {print "$1=>$F[1]"} }' Homo_sapiens.gene_info >entrez.ensembl
4,得到人类entrez ID的go背景
grep '^9606' gene2go |perl -alne '{$hash{$F[1]}.="$F[2] "}END{print "$_\t$hash{$_}" foreach sort keys %hash}' >human.all.go.entrez
5,把人类entrez ID的go背景转换成ensembl的go背景
cat entrez.ensembl human.all.go.entrez |perl -F"=>" -alne  '{$hash{$F[1]}=$F[0];print "$hash{$F[0]}\t$F[1]" if exists $hash{$F[0]}}' >human.all.go.ensembl

在我的群里面共享了所有的代码及帖子内容,欢迎加群201161227,生信菜鸟团!

http://www.bio-info-trainee.com/?p=1

线下交流-生物信息学
同时欢迎下载使用我的手机安卓APP

http://www.cutt.com/app/down/840375

18

查某个基因家族在某物种的具体信息

查某个基因家族在某物种的具体信息

我很伤心,不知道是不是我写的教程还是不够人性化,一个朋友在群里面问如何知道NAC基因家族在拟南芥里面的105个基因信息,我随便给他示范了一下在人类里面如何找,希望他能触类旁通,结果他不会linux,啥生信基础都没有,我只会诱导他简单学习一下,希望他至少明白什么的taxid。所以我给了他我之前写的教程,只希望他告诉我拟南芥的taxid我就帮他把那105个基因找出来。 Continue reading

17

草莓基因组文章解读-并下载原始测序数据

找橡胶测序数据无果

所以我只好找了他们所参考的草莓(strawberry, Fragaria vesca (2n = 2x = 14),a small genome (240 Mb),)的文章,是发表是nature genetics上面的

http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3326587/

可以看到它的SRA索取号。

研读橡胶的基因组文章1087

草莓组装结果: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上面下载这三个数据

研读橡胶的基因组文章1664

 

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

挂在后台自动下载

研读橡胶的基因组文章2877

好了,有了这些数据我们就要进行基因组的一系列分析啦!!!

不过我们可以先看看他们这个研究小组的成果

首先他们建造了一个关于草莓的基因组信息网站

https://strawberry.plantandfood.co.nz/

研读橡胶的基因组文章3091

跟我之前在水科院做鲫鱼鲤鱼的差不多

直接在里面就可以下载他们做好的所有数据,也可以可视化。

 

它的染色体如下,非常简单,就七条染色体

研读橡胶的基因组文章3106

 

http://www.rosaceae.org/species/fragaria/fragaria_vesca/genome_v1.1

我找到了它组装好的草莓基因组地址,用批处理全部下载了

研读橡胶的基因组文章3287

17

研读橡胶的基因组文章-结果没有原始测序数据

研读橡胶的基因组文章

我本科的前两年在海南儋州读书,那时候旁边就是橡胶所,很多同学也在那边做毕业论文什么的,我一直以为那里是全世界的橡胶中心,所有的先进技术都在那里产生,结果,前些天跟一个橡胶所的老师聊天才发现,居然橡胶(Hevea brasiliensis)的基因组已经发表了,可是,跟橡胶所没有半毛钱关系,更搞笑的事情是,堂堂一个基因组文章居然发表在BMC这样的杂志,真不知道是基因组的年代已经过去了还是他们做的实在是太差了,反正我看不过去了,所以研读他们的文章,并且下载数据测试一下。

文章地址如下:http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3575267/

研读橡胶的基因组文章409

可以看到它过于数据的描述都在补充材料1里面,所以我下载了补充材料。

研读橡胶的基因组文章550

可以看到所有的测序数据的描述,45个G的i  llumina的200bp的双端测序,27个G的illumina的200bp的双端测序,约10G左右的长片段(8kb,20kb)罗氏454数据,最后还有一点点solid数据,它这样的测序策略好像是模仿的2011年发布的草莓基因组数据。

 

但是补充材料里面没有列出下载地址,我有点困惑!

按照道理我研读文献的步骤应该没有错,有可能是因为这个文章发表的杂志水平太低,所以不要求他们把测序原始数据上传到NCBI的SRA里面。或者是他们本身觉得文章发的不够档次,不想公布数据,所以先留着自己做精细分析,等发了大文章再公布原始数据。

然后我在NCBI的SRA里面查找了关于橡胶的原始数据,果真没有

研读橡胶的基因组文章727

 

仅有的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上面的

http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3326587/

10

NCBI的taxid简单介绍

 NCBI的taxid简单介绍

物种的信息集合都在它的NCBI的taxid号里面,在NCBI里面关于它的英文介绍地址如下 http://www.ncbi.nlm.nih.gov/guide/taxonomy/ ,NCBI人为的给自然界所有的物种都给了一个编号,这个编号就是taxid,是根据计算机里面树这种数据结构来编码的,其中人类的编号是 9606,7227是果蝇,我们只需要进入这个物种的taxid里面就能看的关于它的一切NCBI存在并且收集好的信息。

NCBI的taxid简单介绍288

 

Continue reading

07

搜索学习其他学者的RNA数据处理流程(包括原始数据、脚本、中间文件)

搜索其他学者的RNA数据处理流程(包括原始数据、脚本、中间文件)

一:原始数据

是谷歌里面无意中搜索到的,是某个物种的RNA数据,不是很大,但是里面有所有的分析流程,非常方便,对原始reads进行了组装,和注释。

http://moana.dnsalias.org/~sgeib/Anth_RNAseq/Run2.1/RawData/

打开网址可以看到raw data的下载链接

QQ截图20150309220349

 

Continue reading

07

阅读文献并下载原始数据知illumina的Chip-seq数据

阅读文献并下载原始数据知illumina的Chip-seq数据

目录

一:阅读文献找到总实验项目

二:在根据实验项目地址找到所有实验数据的下载地址

三:构造脚本并下载

四:用sra-toolkit工具解压

正文

一:阅读文献找到总实验项目

该chip-seq数据其实隶属于一个大的实验项目组,其下载地址如下http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE52964

阅读文献并下载原始数据知illumina的Chip-seq数据482 Continue reading

07

阅读文献并下载原始测序数据之helicos转录组数据

阅读文献并下载原始测序数据之helicos转录组数据

目录

  1. 阅读pdf文献,并找到原始数据搜索关键词。
  2. 根据关键词在NCBI的SRA板块搜索找到其下载地址
  3. 根据下载地址写批处理批量下载所有原始测序数据
  4. 用NCBI提供的工具解压SRR数据,还原成fastq格式reads

正文

一、阅读pdf文献,并找到原始数据搜索关键词

tmp164

可以看到它的下载索引是SRP003040,阅读文献可知其包含4种细胞的6种处理方式的转录组数据

Continue reading