七 19

Bioconductor系列之GenomicAlignments

Posted on by

Bioconductor系列包的安装方法都一样

source("http://bioconductor.org/biocLite.R")biocLite(“GenomicAlignments”)

这个包设计就是用来处理bam格式的比对结果的，具体功能非常多，其实还不如自己写脚本来做这些工作，更方便一点，实在没必要学别人的语法。大家有兴趣可以看GenomicAlignments的主页介绍，三个pdf文档来介绍。

http://www.bioconductor.org/packages/release/bioc/html/GenomicAlignments.html

我们首先读取示例pasillaBamSubset包带有的bam文件

library(pasillaBamSubset)

un1 <- untreated1_chr4()  # single-end reads

library(GenomicAlignments)

reads1 <- readGAlignments(un1)

查看reads1这个结果，可以看到把这个bam文件都读成了一个数据对象GAlignments object，

Bioconductor系列之GenomicAlignments632

readGAlignments这个函数就是读取我们bam文件的，生成的对象可以用多个函数来继续查看信息，比如它的列名都是函数

Seqnames(),strand(),cigar(),qwidth(),start(),end(),width(),njunc() 这些函数对这个GAlignments对象处理后会得到比对的各个情况的向量。

我们也可以读取这个GenomicAlignments包自带的bam文件，而不是pasillaBamSubset包带有的bam文件

> fls <- list.files(system.file("extdata", package="GenomicAlignments"),+                   recursive=TRUE, pattern="*bam$", full=TRUE)

> fls

[1] "C:/Program Files/R/R-3.1.2/library/GenomicAlignments/extdata/sm_treated1.bam"  [2] "C:/Program Files/R/R-3.1.2/library/GenomicAlignments/extdata/sm_untreated1.bam"

reads1 <- readGAlignments(fls[1])

比如我随意截图 cigar(reads1)的结果，就可以看出，大部分都是75M这样的比对情况，可以对这个向量进行统计完全匹配的情况，总共有6077总比对情况。

unique(cigar(reads1))

Bioconductor系列之GenomicAlignments1315

也可以用table格式来查看cigar，可以看到大部分都是M

head(cigarOpTable(cigar(reads1)))

Bioconductor系列之GenomicAlignments1385

接下来我们再读取PE测序数据的比对bam结果，看看分析方法有哪些。

Bioconductor系列之GenomicAlignments1422

U3.galp <- readGAlignmentPairs(untreated3_chr4(), use.names=TRUE, param=param0)head(U3.galp)

也可以分开查看左右两端测序数据的比对情况，跟上面的head是对应关系，上面的

[169,  205] -- [ 326,  362]在下面就被分成左右两个比对情况了。

> head(first(U3.galp))

GAlignments object with 6 alignments and 0 metadata columns:      seqnames strand       cigar    qwidth     start       end     width     njunc         <Rle>  <Rle> <character> <integer> <integer> <integer> <integer> <integer>  [1]     chr4      +         37M        37       169       205        37         0  [2]     chr4      +         37M        37       943       979        37         0  [3]     chr4      +         37M        37       944       980        37         0  [4]     chr4      +         37M        37       946       982        37         0  [5]     chr4      +         37M        37       966      1002        37         0  [6]     chr4      +         37M        37       966      1002        37         0  -------  seqinfo: 8 sequences from an unspecified genome

> head(last(U3.galp))

GAlignments object with 6 alignments and 0 metadata columns:      seqnames strand       cigar    qwidth     start       end     width     njunc         <Rle>  <Rle> <character> <integer> <integer> <integer> <integer> <integer>  [1]     chr4      -         37M        37       326       362        37         0  [2]     chr4      -         37M        37      1086      1122        37         0  [3]     chr4      -         37M        37      1119      1155        37         0  [4]     chr4      -         37M        37       986      1022        37         0  [5]     chr4      -         37M        37      1108      1144        37         0  [6]     chr4      -         37M        37      1114      1150        37         0  -------  seqinfo: 8 sequences from an unspecified genome

用isProperPair可以查看pe测序数据比对情况的，哪些没有配对比对成功

> table(isProperPair(U3.galp))FALSE  TRUE 29518 45828

还可以用Rsamtools包对我们的bam文件进行统计，看下面的代码中fls[1]和un1都是上面得到的bam文件全路径。

> library(Rsamtools)

> quickBamFlagSummary(fls[1], main.groups.only=TRUE)

group |    nb of |    nb of | mean / max                                   of |  records |   unique | records per                              records | in group |   QNAMEs | unique QNAMEAll records........................ A |     1800 |     1800 |    1 / 1  o template has single segment.... S |     1800 |     1800 |    1 / 1  o template has multiple segments. M |        0 |        0 |   NA / NA      - first segment.............. F |        0 |        0 |   NA / NA      - last segment............... L |        0 |        0 |   NA / NA      - other segment.............. O |        0 |        0 |   NA / NA Note that (S, M) is a partitioning of A, and (F, L, O) is a partitioning of M.Indentation reflects this.

> library(Rsamtools)

> quickBamFlagSummary(un1, main.groups.only=TRUE)

group |    nb of |    nb of | mean / max                                   of |  records |   unique | records per                              records | in group |   QNAMEs | unique QNAMEAll records........................ A |   204355 |   190770 | 1.07 / 10  o template has single segment.... S |   204355 |   190770 | 1.07 / 10  o template has multiple segments. M |        0 |        0 |   NA / NA      - first segment.............. F |        0 |        0 |   NA / NA      - last segment............... L |        0 |        0 |   NA / NA      - other segment.............. O |        0 |        0 |   NA / NA Note that (S, M) is a partitioning of A, and (F, L, O) is a partitioning of M.Indentation reflects this.

 

下面讲最后一个综合应用

biocLite("RNAseqData.HNRNPC.bam.chr14")

#这是一个665.5MB的RNA-seq测序数据的比对结果

library(RNAseqData.HNRNPC.bam.chr14)bamfiles <- RNAseqData.HNRNPC.bam.chr14_BAMFILESnames(bamfiles)  # the names of the runslibrary(Rsamtools)quickBamFlagSummary(bamfiles[1], main.groups.only=TRUE)flag1 <- scanBamFlag(isFirstMateRead=TRUE, isSecondMateRead=FALSE,                     isDuplicate=FALSE, isNotPassingQualityControls=FALSE)param1 <- ScanBamParam(flag=flag1, what="seq") #这里是设置readGAlignments的读取参数library(GenomicAlignments)gal1 <- readGAlignments(bamfiles[1], use.names=TRUE, param=param1)

mcols(gal1)$seqoqseq1 <- mcols(gal1)$seqis_on_minus <- as.logical(strand(gal1) == "-")oqseq1[is_on_minus] <- reverseComplement(oqseq1[is_on_minus])

这样就还原得到了原始测序数据啦！

七 19

Bioconductor系列之pasillaBamSubset

Posted on by

这个包主要有bam文件测试数据
> biocLite("pasillaBamSubset")
BioC_mirror: http://bioconductor.orgUsing Bioconductor version 3.0 (BiocInstaller 1.16.5), R version 3.1.2.
Installing package(s) 'pasillaBamSubset'
trying URL 'http://bioconductor.org/packages/3.0/data/experiment/bin/windows/contrib/3.1/pasillaBamSubset_0.3.1.zip'
Content type 'application/zip' length 31514402 bytes (30.1 Mb)
打开pasillaBamSubset包的安装地址就可以看到里面有几个bam文件
Several functions are available for reading BAM files into R:
而且加载包的同时也引入了几个读取bam文件的函数
readGAlignments()
readGAlignmentPairs()
readGAlignmentsList()
scanBam()
Bioconductor系列之pasillaBamSubset448
加载包就可以看到用两个函数得到包自带的数据文件的地址，主要是有很多人不一定把包安装在C盘，所以用函数来定位文件更加安全一点
> library(pasillaBamSubset)
> untreated1_chr4()
[1] "C:/Program Files/R/R-3.1.2/library/pasillaBamSubset/extdata/untreated1_chr4.bam"
> untreated3_chr4()
[1] "C:/Program Files/R/R-3.1.2/library/pasillaBamSubset/extdata/untreated3_chr4.bam"

接下来我们就看看如何读取这些bam文件的
library(pasillaBamSubset)
un1 <- untreated1_chr4() # single-end reads library(GenomicAlignments) reads1 <- readGAlignments(un1) cvg1 <- coverage(reads1) 查看reads1这个结果，可以看到把这个bam文件都读成了一个数据对象GAlignments object， Bioconductor系列之pasillaBamSubset1142
针对着个数据对象有很多操作，其中一个coverage操作是来自于GenomicFeatures
或者GenomicAlignments函数的，可以算出测序覆盖情况。
可以看到这个bam文件里面的比对情况大多几种在4号染色体里面
> cvg1$chr4
integer-Rle of length 1351857 with 122061 runs
Lengths: 891 27 5 12 13 45 5 12 13 ... 5 1 1 3 10
Values : 0 1 2 3 4 5 4 3 2 ... 12 11 10 6
> mean(cvg1$chr4)
[1] 11.33746
> max(cvg1$chr4)[1] 5627
可以看到平均测序深度是11.3X，最大测序深度是5627X

七 18

Bioconductor系列之biomaRt

Posted on by

biomaRt是一个超级网络资源库，里面的信息非常之多，就是网页版的biomaRt的R语言接口。谷歌搜索 the biomart user’s guide filetype:pdf 这个关键词，就看到关于这个包的详细介绍以及例子，我这里简单总结一下它的用法。

这个包一直在更新，函数总是变化，而且需要联网才能使用，基本上等于废物一个，希望大家不要使用它~

包的安装

Continue reading →

七 18

Bioconductor系列之hgu95av2.db

Posted on by

这个包主要都是芯片数据处理方面的应用，本来我是懒得看的，但是我无意中浏览了一下readme，发现它挺适合被作为数据库的典范来学习。
source("http://bioconductor.org/biocLite.R")
biocLite("org.Hs.eg.db")
biocLite("hgu95av2.db")
还是老办法安装了hgu95av2.db之后用ls可以查看这个包里面有36个映射数据，都是把芯片探针ID号转为其它36种主流ID号的映射而已。
library(hgu95av2.db)
> ls("package:hgu95av2.db")
[1] "hgu95av2" "hgu95av2ACCNUM" "hgu95av2ALIAS2PROBE"
[4] "hgu95av2CHR" "hgu95av2CHRLENGTHS" "hgu95av2CHRLOC"
[7] "hgu95av2CHRLOCEND" "hgu95av2.db" "hgu95av2_dbconn"
[10] "hgu95av2_dbfile" "hgu95av2_dbInfo" "hgu95av2_dbschema"
[13] "hgu95av2ENSEMBL" "hgu95av2ENSEMBL2PROBE" "hgu95av2ENTREZID"
[16] "hgu95av2ENZYME" "hgu95av2ENZYME2PROBE" "hgu95av2GENENAME"
[19] "hgu95av2GO" "hgu95av2GO2ALLPROBES" "hgu95av2GO2PROBE"
[22] "hgu95av2MAP" "hgu95av2MAPCOUNTS" "hgu95av2OMIM"
[25] "hgu95av2ORGANISM" "hgu95av2ORGPKG" "hgu95av2PATH"
[28] "hgu95av2PATH2PROBE" "hgu95av2PFAM" "hgu95av2PMID"
[31] "hgu95av2PMID2PROBE" "hgu95av2PROSITE" "hgu95av2REFSEQ"
[34] "hgu95av2SYMBOL" "hgu95av2UNIGENE" "hgu95av2UNIPROT"
但是用这个包自带的函数capture.output(hgu95av2())可以看到这些映射并没有囊括我们标准的hg19版本的2.3万个基因，也就说这个芯片设计的探针只有1.1万个左右。但它根据已有的org.Hs.eg.db来构建了它自己芯片独有的数据包，这样就显得更加正式规范化了，这启发我们研发人员在做一些市场产品的同时也可以把自己的研究成果通主流数据库数据形式结合起来，更易让他人接受。
> capture.output(hgu95av2())
[1] "Quality control information for hgu95av2:"
[2] ""
[3] ""
[4] "This package has the following mappings:"
[5] ""
[6] "hgu95av2ACCNUM has 12625 mapped keys (of 12625 keys)"
[7] "hgu95av2ALIAS2PROBE has 33755 mapped keys (of 103735 keys)"
[8] "hgu95av2CHR has 11540 mapped keys (of 12625 keys)"
[9] "hgu95av2CHRLENGTHS has 93 mapped keys (of 93 keys)"
[10] "hgu95av2CHRLOC has 11474 mapped keys (of 12625 keys)"
[11] "hgu95av2CHRLOCEND has 11474 mapped keys (of 12625 keys)"
[12] "hgu95av2ENSEMBL has 11460 mapped keys (of 12625 keys)"
[13] "hgu95av2ENSEMBL2PROBE has 9814 mapped keys (of 28553 keys)"
[14] "hgu95av2ENTREZID has 11543 mapped keys (of 12625 keys)"
[15] "hgu95av2ENZYME has 2125 mapped keys (of 12625 keys)"
[16] "hgu95av2ENZYME2PROBE has 786 mapped keys (of 975 keys)"
[17] "hgu95av2GENENAME has 11543 mapped keys (of 12625 keys)"
[18] "hgu95av2GO has 11245 mapped keys (of 12625 keys)"
[19] "hgu95av2GO2ALLPROBES has 17214 mapped keys (of 18826 keys)"
[20] "hgu95av2GO2PROBE has 12920 mapped keys (of 14714 keys)"
[21] "hgu95av2MAP has 11519 mapped keys (of 12625 keys)"
[22] "hgu95av2OMIM has 10541 mapped keys (of 12625 keys)"
[23] "hgu95av2PATH has 5374 mapped keys (of 12625 keys)"
[24] "hgu95av2PATH2PROBE has 228 mapped keys (of 229 keys)"
[25] "hgu95av2PMID has 11529 mapped keys (of 12625 keys)"
[26] "hgu95av2PMID2PROBE has 372382 mapped keys (of 432400 keys)"
[27] "hgu95av2REFSEQ has 11506 mapped keys (of 12625 keys)"
[28] "hgu95av2SYMBOL has 11543 mapped keys (of 12625 keys)"
[29] "hgu95av2UNIGENE has 11533 mapped keys (of 12625 keys)"
[30] "hgu95av2UNIPROT has 11315 mapped keys (of 12625 keys)"
[31] ""
[32] ""
[33] "Additional Information about this package:"
[34] ""
[35] "DB schema: HUMANCHIP_DB"
[36] "DB schema version: 2.1"
[37] "Organism: Homo sapiens"
[38] "Date for NCBI data: 2014-Sep19"
[39] "Date for GO data: 20140913"
[40] "Date for KEGG data: 2011-Mar15"
[41] "Date for Golden Path data: 2010-Mar22"
[42] "Date for Ensembl data: 2014-Aug6"
这个hgu95av2.db所加载的36个包都是一种特殊的对象，但是可以把它当做list来操作，是一种类似于hash的对应表格，其中keys是独特的，但是value可以有多个。
既然是类似于list，那我就简单讲几个小技巧来操作这些数据对象。所有的操作都要用as.list()函数来把数据展现出来
> as.list(hgu95av2ENZYME[1])
$`1000_at`
[1] "2.7.11.24"
可以看到这样就提取出来了hgu95av2ENZYME的第一个元素，key是`1000_at`，它所映射的value是 "2.7.11.24"这个酶。
同理对list取元素的三个操作在这里都可以用
> as.list(hgu95av2ENZYME['1000_at'])
$`1000_at`
[1] "2.7.11.24"
> as.list(hgu95av2ENZYME$'1000_at')
[[1]]
[1] "2.7.11.24"
然而这只是list的操作，我们还有一个函数专门对这个对象来取元素，那就是get函数，取多个元素用mget，类似于as.list(hgu95av2ENZYME[1:10])
用函数mget(probes_id,hgu95av2ENZYME)就可以根据自己的probes_id这个向量来选取数据了，当然也要用as.list()来展示出来。
值得一提的是这里有个非常重要的函数是revmap，可以把我们的key和value调换位置。
因为我们的数据对象里面的对应关系不是一对一，而是一对多，例如一个基因往往有多个go通路信息，例如这个就有十几个
as.list(hgu95av2GO$'1000_at')
$`GO:0000165`
$`GO:0000165`$GOID
[1] "GO:0000165"

$`GO:0000165`$Evidence
[1] "TAS"

$`GO:0000165`$Ontology
[1] "BP"

$`GO:0000186`
$`GO:0000186`$GOID
[1] "GO:0000186"

$`GO:0000186`$Evidence
[1] "TAS"

$`GO:0000186`$Ontology
[1] "BP"
一层层的数据结构非常清晰，但是数据太多，所以它给了一个toTable函数来格式化，就是把一对多的list压缩成表格。
> toTable(hgu95av2GO[1])
probe_id go_id Evidence Ontology
1 1000_at GO:0000165 TAS BP
2 1000_at GO:0000186 TAS BP
3 1000_at GO:0000187 TAS BP
4 1000_at GO:0002224 TAS BP
5 1000_at GO:0002755 TAS BP
6 1000_at GO:0002756 TAS BP
7 1000_at GO:0006360 TAS BP
------------------------------------------------------------
81 1000_at GO:0004707 NAS MF
82 1000_at GO:0001784 IEA MF
83 1000_at GO:0005524 IEA MF
这样就非常方便的查看啦。
当然对这个数据对象还有很多实用的操作。length(),Llength(),Rlength(),keys(),Lkeys(),Rkeys()等等，还有count.mappedkeys(),mappedLkeys()，还有个比较特殊的isNA来判断是否有些keys探针没有对应任何信息。
而且以上这些函数不仅可以用来获取数据对象的信息，还可以修改这个数据对象。
Lkeys(hgu95av2CHR) 可以查看这个对象有12625个探针。
而> Rkeys(hgu95av2CHR)
[1] "19" "12" "8" "14" "3" "2" "17" "16" "9" "X" "6" "1" "7" "10" "11"
[16] "22" "5" "18" "15" "Y" "20" "21" "4" "13" "MT" "Un"
可以看到这些探针分布在不同的染色体上面，如果我这时候给
x=hgu95av2CHR
Rkeys(x) = c("1","2")
toTable(x)可以看到这时候x只剩下1946个探针了，就是1,2号染色体上面的基因了。
然后可以一个个看这些函数的用法，其实就是SQL的增删改查的基本操作而已。
> length(x)
[1] 12625
> Llength(x)
[1] 12625
> Rlength(x)
[1] 2
> head(keys(x))
[1] "1000_at" "1001_at" "1002_f_at" "1003_s_at" "1004_at" "1005_at"
> head(Lkeys(x))
[1] "1000_at" "1001_at" "1002_f_at" "1003_s_at" "1004_at" "1005_at"
> head(Rkeys(x))
[1] "1" "2"
> count
count.fields count.mappedLkeys countOverlaps counts<- count.links count.mappedRkeys countQueryHits countSubjectHits count.mappedkeys countMatches counts > count.mappedkeys(x)
[1] 1946
> count.mappedLkeys(x)
[1] 1946
> count.mappedRkeys(x)
[1] 2

七 17

perl操作pdf文档

Posted on by

大家看看就好，这个模块写的不怎么样，而且有高手已经写了一个pdftoolkit就是完全用这个模块实现了大部分pdf文档的操作

PDF::API2模块使用笔记

一：简单使用方法

use PDF::API2;

# Create a blank PDF file $pdf = PDF::API2->new();
# Open an existing PDF file $pdf = PDF::API2->open('some.pdf');
# Add a blank page $page = $pdf->page();
# Retrieve an existing page $page = $pdf->openpage($page_number);
# Set the page size $page->mediabox('Letter');
# Add a built-in font to the PDF $font = $pdf->corefont('Helvetica-Bold');
# Add an external TTF font to the PDF $font = $pdf->ttfont('/path/to/font.ttf');
 

# Add some text to the page

 $text = $page->text();

$text->font($font, 20);

$text->translate(200, 700);

$text->text('Hello World!');
# Save the PDF $pdf->saveas('/path/to/new.pdf');

 

实例：

use PDF::API2;

$pdf=PDF::API2->new;

$pdf->mediabox('A4');

$ft=$pdf->cjkfont('Song');

$page = $pdf->page;

$gfx=$page->gfx;

$gfx->textlabel(50,750,$ft,20,"\x{Cool44}\x{4EA7}"); # 资产二字

$pdf->saveas('Song_Test.pdf');

 

二：主要对象及方法

1、pdf对象可以创造，可以打开，可以保存，可以更新，还有一堆参数可以设置

$pdf->preferences(%options)还可以设置一些浏览参数，不过本来pdf阅读器可以设置，没必要在这里花时间。

这个可以当做是个人创建pdf的保密信息，也许有一点用吧。

还可以可以设置页脚$pdf->pageLabel($index, $options

2、Page对象，可以新建，可以打开，可以保存（需要指定保存的位置）

$page = $pdf->page()

$page = $pdf->page($page_number)

$page = $pdf->openpage($page_number);

还可以更新旧的pdf，这样可以循环获取pdf页面不停的累积到一个新的pdf

$page = $pdf->import_page($source_pdf, $source_page_number, $target_page_number)

$pdf = PDF::API2->new();

$old = PDF::API2->open('our/old.pdf');   # Add page 2 from the old PDF as page 1 of the new PDF

$page = $pdf->import_page($old, 2);

$pdf->saveas('our/new.pdf');If $source_page_number is 0 or -1, it will return the last page in the document.

$count = $pdf->pages()Returns the number of pages in the document.

这样就可以写一个简单程序把我们的pdf文件合并

use PDF::API2;

my $new = PDF::API2->new;

foreach my $filename (@ARGV) {   my $pdf = PDF::API2->open($filename);   $new->importpage($pdf, $_) foreach 1 .. $pdf->pages;}$new->saveas('new.pdf'); $pdf->mediabox($name)

可以指定A4，A3，A5等等$pdf->mediabox($w, $h)可以指定宽度和高度$pdf->mediabox($llx, $lly, $urx, $ury)

3，还可以随意画点线面及表格，太复杂了就不看了

 

七 17

转-windows快捷键，让你的办公效率提升一个档次

Posted on by
  1. gpedit.msc-----组策略

2. sndrec32-------录音机

3. Nslookup-------IP地址侦测器

4. explorer-------打开资源管理器

5. logoff---------注销命令

6. tsshutdn-------60秒倒计时关机命令

7. lusrmgr.msc----本机用户和组

8. services.msc---本地服务设置

9. oobe/msoobe /a----检查XP是否激活

10. notepad--------打开记事本

11. cleanmgr-------垃圾整理

12. net start messenger----开始信使服务

13. compmgmt.msc---计算机管理

15. conf-----------启动netmeeting

16. dvdplay--------DVD播放器

17. charmap--------启动字符映射表

18. diskmgmt.msc---磁盘管理实用程序

　19. calc-----------启动计算器

20. dfrg.msc-------磁盘碎片整理程序

21. chkdsk.exe-----Chkdsk磁盘检查

22. devmgmt.msc--- 设备管理器

23. regsvr32 /u *.dll----停止dll文件运行

24. drwtsn32------ 系统医生

25. rononce -p ----15秒关机

26. dxdiag---------检查DirectX信息

27. regedt32-------注册表编辑器

28. Msconfig.exe---系统配置实用程序

29. rsop.msc-------组策略结果集

30. mem.exe--------显示内存使用情况

31. regedit.exe----注册表

32. winchat--------XP自带局域网聊天

33. progman--------程序管理器

34. winmsd---------系统信息

　　43. write----------写字板

44. winmsd---------系统信息

46. winchat--------XP自带局域网聊天

48. Msconfig.exe---系统配置实用程序

49. mplayer2-------简易widnows media player

50. mspaint--------画图板

51. mstsc----------远程桌面连接

52. mplayer2-------媒体播放机

53. magnify--------放大镜实用程序

　54. mmc------------打开控制台

55. mobsync--------同步命令

56. dxdiag---------检查DirectX信息

57. drwtsn32------ 系统医生

58. devmgmt.msc--- 设备管理器

59. dfrg.msc-------磁盘碎片整理程序

60. diskmgmt.msc---磁盘管理实用程序

61. dcomcnfg-------打开系统组件服务

62. ddeshare-------打开DDE共享设置

65. net start messenger----开始信使服务

67. nslookup-------网络管理的工具向导

68. ntbackup-------系统备份和还原

69. narrator-------屏幕“讲述人”

70. ntmsmgr.msc----移动存储管理器

71. ntmsoprq.msc---移动存储管理员操作请求

72. netstat -an----(TC)命令检查接口

73. syncapp--------创建一个公文包

　　74. sysedit--------系统配置编辑器

75. sigverif-------文件签名验证程序

76. sndrec32-------录音机

77. shrpubw--------创建共享文件夹

78. secpol.msc-----本地安全策略

　80. services.msc---本地服务设置

81. Sndvol32-------音量控制程序

82. sfc.exe--------系统文件检查器

83. sfc /scannow---windows文件保护

84. tsshutdn-------60秒倒计时关机命令

85. tourstart------xp简介(安装完成后出现的漫游xp程序)

86. taskmgr--------任务管理器

　87. eventvwr-------事件查看器

88. eudcedit-------造字程序

　92. progman--------程序管理器

94. rsop.msc-------组策略结果集

95. regedt32-------注册表编辑器

96. rononce -p ----15秒关机

99. cmd.exe--------CMD命令提示符

100. chkdsk.exe-----Chkdsk磁盘检查

101. certmgr.msc----证书管理实用程序

　102. calc-----------启动计算器

103. charmap--------启动字符映射表

104. cliconfg-------SQL SERVER 客户端网络实用程序

105. Clipbrd--------剪贴板查看器

106. conf-----------启动netmeeting

107. compmgmt.msc---计算机管理

108. cleanmgr-------垃圾整理

109. ciadv.msc------索引服务程序

110. osk------------打开屏幕键盘

　113. lusrmgr.msc----本机用户和组

114. logoff---------注销命令

115. fsmgmt.msc-----共享文件夹管理器

116. utilman--------辅助工具管理器

117. iexpress-------木马捆绑工具

　打开服务管理器的是services.msc

如果要用cmd直接启用已知服务名的服务如下：

net start [服务名] 启动一个服务

net stop [服务名] 停用一个服务

 

七 17

读书笔记-核酸&蛋白序列分析

Posted on by

核酸序列分析：
一．DNA基本序列分析：bioedit、DNAman、DNAstar
a)    组成成分分析
b)    序列转换分析
c)    酶切位点分析（REBSASE，NEBCutter2）
二．DNA序列特征分析
a)    开放阅读框分析
b)    启动子和转录因子结合位点（数据库EPD,TRANSFAC,DBTSS,TRRD）（软件Promoter Scan，TfBlast，TESS）
c)    CpG岛识别（在线工具：EMBOSS,CpG Island searcher，CpG cluster）
三．重复序列分析
a)    常用数据库（RepBase，RepeatMasker，LINE-1，STR）
四．基因识别，
a)    同源序列比对
b)    组成统计学特征预测（UTR,EXON,INTRON,）
c)    GENESCAN，GRAIL，geneMarkS，Glimmer
五．mRNA可变剪切分析
a)    常用数据库（ASTD,ASD,ASAP）
b)    在线工具（ASPicDB,ESEfinder,RESCUE-ESE）
六．miRNA与靶基因预测
a)    预测方法（同源片段搜索，比较基因组学预测，序列结构特征打分，靶标预测，机器学习）
b)    数据库资源（miRNA，miRBase，TarBase，miRGen，MIRSCAN，MiPred，miRFinder，miRanda，Targetscan，PicTar）

蛋白序列分析
一．蛋白基本序列分析
a)    氨基酸组成（ExPASy）
b)    电荷性质，疏水性（ProtScale）
c)    理化性质（ProtParam）
二．序列特征信息
a)    跨膜区预测（TMpred）
b)    信号肽分析（signalP）
c)    卷曲螺旋分析（COILS）
三．功能信息分析
a)    PROSITE蛋白质功能数据库
b)    结构域和功能位点分析InterProScan
c)    基于蛋白同源性的功能分析（blastp等）
四．蛋白质结构分析
a)    二级结构预测-（Chou-Fasman，GOR，PHD，NNSSP以及多元预测方法）
b)    三级结构预测-（同源建模。重头预测）

七 17

用DESeq进行差异分析的源代码

Posted on by

要保证当前文件夹下面有了742KO1.count等4个文件，就是用htseq对比对的bam文件进行处理后的输出文件

library(DESeq)
#加载数据
K1=read.table("742KO1.count",row.names=1)
K2=read.table("743KO2.count",row.names=1)
W1=read.table("740WT1.count",row.names=1)
W2=read.table("741WT2.count",row.names=1)
#列名
data=cbind(K1,K2,W1,W2)
#如果是htseq的结果，则删除data最后四行
n=nrow(data)
data=data

[c language="(-n+4:-n),"][/c]

#如果是bedtools的结果，取出统计个数列和行名
kk1=cbind(K1$V5)
rownames(kk1)=rownames(K1)
K1=kk1

#差异分析
colnames(data)=c("K1","K2","W1","W2")
type=rep(c("K","W"),c(2,2))
de=newCountDataSet(data,type)
de=estimateSizeFactors(de)
de=estimateDispersions(de)
res=nbinomTest(de,"K","W")

#res就是我们的表达量检验结果

到这里，理论上差异基因的分析已经结束啦！后面只是关于R的bioconductor包的一些简单结合使用而已

library(org.Mm.eg.db)

tmp=select(org.Mm.eg.db, keys=res$id, columns=c("ENTREZID","SYMBOL"), keytype="ENSEMBL")

#合并res和tmp
res=merge(tmp,res,by.x="ENSEMBL",by.y="id",all=TRUE)

#go
tmp=select(org.Mm.eg.db, keys=res$ENSEMBL, columns="GO", keytype="ENSEMBL")
ensembl_go=unlist(tapply(tmp[,2],as.factor(tmp[,1]),function(x) paste(x,collapse ="|"),simplify =F))

#为res加入go注释，
res$go=ensembl_go[res$ENSEMBL]#为res加入一列go

#写入all——data
all_res=res
write.csv(res,file="all_data.csv",row.names =F)

uniq=na.omit(res)#删除无效基因
sort_uniq=uniq[order(uniq$padj),]#按照矫正p值排序

#写入排序后的all_data
write.csv(res,file="all_data.csv",row.names =F)

#标记上下调基因
sort_uniq$up_down=ifelse(sort_uniq$baseMeanA>sort_uniq$baseMeanB,"up","down")
#交换上下调基因列位置
final_res=sort_uniq[,c(12,1:11)]
#写出最后数据
write.csv(final_res,file="final_annotation_gene_bedtools_sort_uniq.csv",row.names =F)

#然后挑选出padj值小于0.05的数据来做富集
tmp=select(org.Mm.eg.db, keys=sort_uniq[sort_uniq$padj<0.05,1], columns="ENTREZID", keytype="ENSEMBL")
diff_ENTREZID=tmp$ENTREZID
require(DOSE)
require(clusterProfiler)
diff_ENTREZID=na.omit(diff_ENTREZID)
ego <- enrichGO(gene=diff_ENTREZID,organism="mouse",ont="CC",pvalueCutoff=0.05,readable=TRUE)
ekk <- enrichKEGG(gene=diff_ENTREZID,organism="mouse",pvalueCutoff=0.05,readable=TRUE)
write.csv(summary(ekk),"KEGG-enrich.csv",row.names =F)
write.csv(summary(ego),"GO-enrich.csv",row.names =F)

 

七 17

转载-VCF格式详解

Posted on by

CHROM(chromosome):染色体

POS - position:参考基因组variant碱基位置，如果是INDEL（插入缺失），位置是INDEL的第一个碱基位置

ID - identifier: variant的ID。比如在dbSNP中有该SNP的id，则会在此行给出；若没有，则用’.'表示其为一个novel variant。

REF - reference base(s):参考碱基，染色体上面的碱基，必须是ATCGN中的一个，N表示不确定碱基

ALT - alternate base(s):与参考序列比较发生突变的碱基

QUAL - quality: Phred格式(Phred_scaled)的质量值，表 示在该位点存在variant的可能性；该值越高，则variant的可能性越大；计算方法：Phred值 = -10 * log (1-p) p为variant存在的概率; 通过计算公式可以看出值为10的表示错误概率为0.1，该位点为variant的概率为90%。

FILTER - _filter status: 使用上一个QUAL值来进行过滤的话，是不够的。GATK能使用其它的方法来进行过滤，过滤结果中通过则该值为”PASS”;若variant不可靠，则该项不为”PASS”或”.”。

INFO - additional information:  这一行是variant的详细信息，具体如下：

DP-read depth：样本在这个位置的reads覆盖度。是一些reads被过滤掉后的覆盖度。          DP4:高质量测序碱基，位于REF或者ALT前后

MQ：表示覆盖序列质量的均方值RMS Mapping Quality

FQ：phred值关于所有样本相似的可能性

AF1： AF(Allele Frequency) 表示Allele的频率，AF1为第一个ALT allele 发生频率的可能性评估

AC1：AC表示Allele（等位基因）的数目,AC1为对第一个ALT allele count的最大可能性评估

AN：AN(Allele Number) 表示Allele的总数目

IS：插入缺失或部分插入缺失的reads允许的最大数量

AC：AC(Allele Count) 表示该Allele的数目

G3：ML 评估基因型出现的频率

HWE：chi^2基于HWE的测试p值和G3

CLR：在受到或者不受限制的情况下基因型出现可能性log值

UGT：最可能不受限制的三种基因型结构

CGT：最可能受限制三种基因型的结构

PV4：四种P值得误差，分别是（strand、baseQ、mapQ、tail distance bias）

INDEL：表示该位置的变异是插入缺失

PC2：非参考等位基因的phred（变异的可能性）值在两个分组中大小不同

PCHI2：后加权chi^2，根据p值来测试两组样本之间的联系

QCHI2：Phred scaled PCHI2.

PR：置换产生的一个较小的PCHI2

QBD：Quality by Depth，测序深度对质量的影响

RPB：序列的误差位置（Read Position Bias）

MDV：样本中高质量非参考序列的最大数目

VDB：Variant Distance Bias，RNA序列中过滤人工拼接序列的变异误差范围

GT：样品的基因型（genotype）。两个数字中间用’/'分 开，这两个数字表示双倍体的sample的基因型。0 表示样品中有ref的allele； 1 表示样品中variant的allele； 2表示有第二个variant的allele。因此： 0/0 表示sample中该位点为纯合的，和ref一致； 0/1 表示sample中该位点为杂合的，有ref和variant两个基因型； 1/1 表示sample中该位点为纯合的，和variant一致。

GQ：基因型的质量值(Genotype Quality)。Phred格式(Phred_scaled)的质量值，表示在该位点该基因型存在的可能性；该值越高，则Genotype的可能性越 大；计算方法：Phred值 = -10 * log (1-p) p为基因型存在的概率。

GL：三种基因型（RR RA AA）出现的可能性，R表示参考碱基，A表示变异碱基

DV：高质量的非参考碱基

SP：phred的p值误差线

PL：指定的三种基因型的质量值(provieds the likelihoods of the given genotypes)。这三种指定的基因型为(0/0,0/1,1/1)，这三种基因型的概率总和为1。和之前不一致，该值越大，表明为该种基因型的可能 性越小。 Phred值 = -10 * log (p) p为基因型存在的概率。

FORMAT 和 BC1-1-base.sorted.bam：这两行合起来提供了’ BC1-1-base′这个sample的基因型的信息。’ BC1-1-base′代表这该名称的样品，是由BAM文件中的@RG下的 SM 标签决定的。

七 17

转-基因突变种类大全

Posted on by

突变（Mutation, 即基因突变）:在生物学上的含义，是指细胞中的遗传基因（通常指存在于细胞核中的脱氧核糖核酸）发生的改变。它包括单个碱基改变所引起的点突变，或多个碱基的缺失、重复和插入。原因可以是细胞分裂时遗传基因的复制发生错误、或受化学物质、辐射或病毒的影响。

以功能分类:

失去功能的突变Loss-of-function mutations

：失去功能的突变是指发生的突变会造成基因完全地失去活性，原因可分成两类。一类是由于基因被删除或是调控基因表现的过程受到影响让基因不表现，另一种则是由于基因本身受到影响，使得基因的产物蛋白质失去功能。又称剔除突变（null mutations）或是敲除突变（knockout mutations）。

次形态突变Time form mutation：此种突变会使基因的表现或是基因产物的活性减弱，但不会消失。

超形态突变hypermorphic mutations：此种突变与次形态突变相反，会使基因的表现加强。

获得功能的突变gain-of-function mutation：获得功能的突变是指发生的突变让原本应该是不表现的基因产生活性，进而影响细胞功能，这样的突变多半需要染色体程度的突变较有可能产生，而最常发生获得功能的突变就是癌细胞。

以突变机理分类:

  1. 点突变point mutation:DNA序列中涉及单个核苷酸或碱基的变化称为点突变。 通常有两种情况：一是一种碱基或核苷酸被另一种碱基或核苷酸所替换；二是一个碱基的插入或缺失。

                   （1）沉默突变silent mutation

当点突变发生在基因及其调控序列之外，或使基因序列内一种密码子变成编码同一种氨基酸的另一种同义密码子时，不会改变生物个体的基因产物，因而不引起性状变异。不引起生物性状变异的突变称为沉默突变。

                   （2）错义突变missense mutation

指由于某个碱基对的改变，使编码一种氨基酸的密码子变成编码另外一种氨基酸的密码子，结果使构成蛋白质的数百上千个氨基酸中有一个氨基酸发生变化。（实例：镰刀形细胞贫血症）

                   （3）移码突变frameshift mutation

指在DNA链上，有时一个或几个非3的整数倍的碱基的插入或缺失，往往产生比碱基替换突变更严重的后果。 这种插入或缺失突变会造成阅读框的改变，翻译过程中其下游的三联密码子都被错读，产生完全错误的肽链或肽链合成提前终止。这种插入或缺失突变又称为移码突变。

                   （4）无义突变nonsense mutation

是指当点突变使一个编码氨基酸的密码子变成终止子时，则蛋白质合成进行到该突变位点时会提前终止，结果产生一个较短的多肽链或较小的蛋白质。

  1. 大突变

大突变是可能涉及整个基因以至多个基因的一长段DNA序列的改变，大突变常常导致染色体畸变。

（1）缺失:指DNA分子丢失一段碱基序列。（染色体缺失）(Deletion)

（2）插入:指DNA分子的正常序列中插入一段DNA序列。(Insertion②)

（3）重排:重排包括某段DNA序列的重复（duplication），倒位（inversion），易位（translocation）等。

 

七 17

转-R语言内存管理

Posted on by

R中的对象（比如矩阵）在内存中存于两种不同的地方:

第一种是堆内存（heap），其基本单元是“Vcells”，每个大小为8字节，新来一个对象就会申请一块空间，把值全部存在这里，和C里面的堆内存很像;

第二种是地址对（cons cells），主要用来存储地址信息，最小单元一般在32位系统中是28字节、64位系统中是56字节。

1、ls()来查看当前所有对象名，对于每一个对象，可以通过object.size(x)来查看其占用内存的大小。

如果是因为当前对象占用内存过多，那么可以通过处理对象来获取更大的可用内存。一个很有用的方法是改变对象的存储模式，通过storage.mode（x）可以看到某个对象的存储模式，比如某个矩阵默认就是“double”的，如果这个矩阵的数值都是整数甚至0-1，完全没必要使用double来占用空间，可以使用storage.mode(x） <- "integer"将其改为整数型，可以看到该对象的大小会变为原来的一半。

2、object.size()看每个变量占多大内存。

3、memory.size()查看现在的work space的内存使用

4、memory.limit()查看系统规定的内存使用上限。如果现在的内存上限不够用，可以通过memory.limit(newLimit)更改到一个新的上限。注意，在32位的R中，封顶上限为4G，无法在一个程序上使用超过4G （数位上限）。这种时候，可以考虑使用64位的版本。

 

对于一些很大的但无用的中间变量，养成清理的习惯：

可以使用rm(object)删除变量，但是记住，rm后记得使用gc()做Garbage collection，否则内存是不会自动释放的，相当于你没做rm.

七 16

模拟测序lambda_virus基因组

Posted on by
lambda_virus基因组文件是bowtie软件自带的测试数据，共48502个bp，首先我用脚本模拟出它的全打断文件！

perl -alne '{next if /^>/;$a.=$_;}END{$len=length $a;print substr($a.$a,$_,120) foreach 0..$len}' lambda_virus.fa >lamb_virus.120bp

长度均为120bp的片段。

我测序的策略是CTAG碱基重复30次，共加入120个碱基。

对每个120bp片段来说，如果遇到互补碱基就加上，直到120个碱基加完，这样如果比较巧合的话，会有部分碱基能全部加满120bp的，但是如果每个120bp片段的ATCG分布均匀，那么就都应该30bp碱基能被加上。

image001

[perl]
while (<>) {

$seq=$_;$sum=0;

foreach $i (0..120){

$str=substr($seq,$i,2);

if ($str eq "GG"| $str eq "CC"| $str eq "AA"| $str eq "TT"){$sum+=4;}

elsif ($str eq "GT"| $str eq "CG"| $str eq "AC"| $str eq "TA"){$sum+=3;}

elsif ($str eq "GA"|$str eq "CT"| $str eq "AG"| $str eq "TC"){$sum+=2;}

else{$sum+=1;};

#print "$sum\n";

if ($sum>120){print "$i\n";last;}

}

}

[/perl]

perl length.pl lambda_virus.120bp >length.txt

得到结果如下：

 

Length 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51
Count 2 19 34 110 204 432 878 1495 2237 3202 4343 5179 5697 5429 4865 4214
Length 52 53 54 55 56 57 58 59 60 61 62 63 64
Count 3249 2499 1735 1090 657 396 228 141 90 48 18 9 3

右表可以看出，大部分测序得到碱基长度都集中在46bp到51bp之间

画出箱线图如下

image003

画出条形图如下：

image005

 

然后我模拟了一个6000bp的基因组，做同样的模拟测序看看评价测序长度分布情况：

Length 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58
Count 9 22 96 207 322 382 479 671 770 714 706 546 424 232 182 100 52 30 14 21
Length 59 60 61                                  
Count 15 5 2                                  

可以看出这次的测序片段集中在45到51bp

七 15

perl操作excel表格

Posted on by

perl这个语言已经过时很久了，所以它的模块支持性能不是很好，暂时我只看到了对excel2003格式的表格的读取写入操作。

以下是我参考Spreadsheet::ParseExcel这个模块写的一个把excel表格转为csv的小程序，大家也可以自己搜索该模块的说明文档，这样学的更快一点！

[perl]

#!/usr/bin/perl -w
# For each tab (worksheet) in a file (workbook),
# spit out columns separated by ",",
# and rows separated by c/r.

use Spreadsheet::ParseExcel;
use strict;
use utf8;
use Encode::Locale qw($ENCODING_LOCALE_FS);
use Encode;

my $filename ="test.xls";#输入需要解析的excel文件名，必须是03版本的
my $e = new Spreadsheet::ParseExcel; #新建一个excel表格操作器
my $eBook = $e->Parse($filename);    #用表格操作器来解析我们的文件
my $sheets = $eBook->{SheetCount};   #得到该文件中sheet总数
my ($eSheet, $sheetName);

foreach my $sheet (0 .. $sheets - 1) {
$eSheet = $eBook->{Worksheet}[$sheet];
$sheetName = $eSheet->{Name};
my $f1 = encode(locale_fs => $sheetName); #每次操作中文我都很纠结，还得各种转码
open FH_out ,">$f1.csv" or die "error open ";
next unless (exists ($eSheet->{MaxRow}) and (exists ($eSheet->{MaxCol})));
foreach my $row ($eSheet->{MinRow} .. $eSheet->{MaxRow}) {
foreach my $column ($eSheet->{MinCol} .. $eSheet->{MaxCol}) {
if (defined $eSheet->{Cells}[$row][$column])
{
print FH_out $eSheet->{Cells}[$row][$column]->Value . ",";
} else {
print FH_out ",";
}
}
print FH_out "\n";
}
close FH_out;

}
exit;

[/perl]

 

七 13

使用Seq2HLA进行HLA分型

Posted on by

基于高通量测序数据进行HLA分型的软件挺多的，比较老的有三个，作者分别是Boegel et al.Kim et al.Major et al.，然后他们都被OptiType这个软件的作者被批评了，我这里先介绍Kim et al的seq2HLA使用方法，以下是它的一些链接。

功能概述

seq2HLA is a computational tool to determine Human Leukocyte Antigen (HLA) directly from existing and future short RNA-Seq reads. It takes standard RNA-Seq sequence reads in fastq format as input, uses a bowtie index comprising known HLA alleles and outputs the most likely HLA class I and class II types, a p-value for each call, and the expression of each class.

软件简介

Type of tool     Program

Nature of tool          Standalone

Operating system   Unix/Linux, Mac OS X

Language        Python, R

Article     (Boegel et al., 2013) HLA typing from RNA-Seq sequence reads. Genome medicine.

PubMed http://www.ncbi.nlm.nih.gov/pubmed/23259685

URL          https://bitbucket.org/sebastian_boegel/seq2hla

源代码，下载并安装

https://bitbucket.org/sebastian_boegel/seq2hla/src

http://tron-mainz.de/tron-facilities/computational-medicine/seq2hla/

第一版是这样的

image001

第二版是这样的

image002

只有第二版才支持gz压缩包格式的fastq，而且不需要指定length了

其中reference文件夹下面的是发布这个软件的团体已经制备好来的HLA库文件

image003

下载即可使用，前提是你的系统其它环境都OK

用法：

python seq2HLA.py -1 <readfile1> -2 <readfile2> -r "<runname>" [-p <int>]* [-3 <int>]**

image004

很简单，-1和-2指定我们的双端测序数据即可，可以是压缩包格式的（自动调用gzip），-r的输出目录，会输出7个文件，需要一个个解读，-p指定线程数给bowtie用的，-3是指定需要trim几个低质量碱基。

但是运行这个软件的要求非常多，需要安装好python和R，而且还有版本限制，需要安装好biopython而且还必须是双端测序，而且当前文件夹下面的reference文件夹下面必须有参考基因组的bowtie索引，而且系统必须安装好了bowtie，还需要在快捷方式里面！

我这里用的是第二版的

image006

所以，我用的也是第二版改进的命令。非常好用,我这里用的是一个外显子测序数据，是hiseq2500测的PE100

python seq2HLA.py -1 ../../6-exon/PC3-1.read1_Clean.fastq.gz -2 ../../6-exon/PC3-1.read2_Clean.fastq.gz -r PC3

貌似输出文件太多了一点

#Output:#The results are output to stdout and to textfiles. Most important are:

#i) <prefix>-ClassI.HLAgenotype2digits => 2 digit result of Class I

#ii) <prefix>-ClassII.HLAgenotype2digits => 2 digit result of Class II

#iii) <prefix>-ClassI.HLAgenotype4digits => 4 digit result of Class I

#iv) <prefix>-ClassII.HLAgenotype4digits => 4 digit result of Class II

#v) <prefix>.ambiguity => reports typing ambuigities (more than one solution for an allele possible)

#vi) <prefix>-ClassI.expression => expression of Class I alleles

#vii) <prefix>-ClassII.expression => expression of Class II alleles

根据文献，我简单看了一下，文件的确好复杂，不过我们只需要看输出日志即可

-----------2 digit typing results-------------

#Locus Allele 1       Confidence     Allele 2       Confidence

A       A*68   7.287148e-05   A*24   0.03680272

B       B*52   0.1717737       B*53   0.3952319

C       C*12   0.03009331     hoz("C*14")     0.6783964

Calculation of locus-specific expression ...

BC1-1/BC1-1-ClassI.bowtielog

A: 7.93 RPKM

C: 9.75 RPKM

B: 8.35 RPKM

The digital haplotype is written into BC1-1/BC1-1-ClassI.digitalhaplotype3

-----------4 digit typing results-------------

#Locus Allele 1       Confidence     Allele 2       Confidence

!A     A*68:01 7.287148e-05   A*24:02 0.03680272

!B     B*52:01 0.1717737       B*53:01'       0.6542288

!C     C*12:02 0.03371717     C*12:02 0.6783964

上面的HLA的class I的数据结果

接下来是class II的数据结果，是不是很简单呀！

-----------2 digit typing results-------------

#Locus Allele 1       Confidence     Allele 2       Confidence

DQA     DQA1*01 0.1511134       DQA1*02 0

DQB     DQB1*02 0.02321615     DQB1*05 0.42202

DRB     DRB1*15 2.595144e-05   DRB1*07 0.321219

Calculation of locus-specific expression ...

BC1-1/BC1-1-ClassII.bowtielog

DQB1: 4.47 RPKM

DRB1: 5.59 RPKM

DQA1: 0.44 RPKM

-----------4 digit typing results-------------

#Locus Allele 1       Confidence     Allele 2       Confidence

!DQA   DQA1*01:02'     0.1511134       DQA1*02:01     0.0

!DQB   DQB1*02:01'     0.02321615     DQB1*05:01     0.42202

!DRB   DRB1*15:02'     2.595144e-05   DRB1*07:01     0.321219

七 06

GATK使用注意事项

Posted on by

GATK这个软件在做snp-calling的时候使用率非常高，因为之前一直是简单粗略的看看snp情况而已，所以没有具体研究它。

这些天做一些外显子项目以找snp为重点，所以想了想还是用起它，报错非常多，调试了好久才成功。

所以记录一些注意事项!

GATK软件本身是受版权保护的，所以需要申请才能下载使用，大家自己去broad institute申请即可。

下载软件就可以直接使用，java软件不需要安装，但是需要你的机器上面有java，当然软件只是个开始，重点是你还得下载很多配套数据，https://software.broadinstitute.org/gatk/download/bundle（ps:这个链接可能会失效，下面的文件，请自己谷歌找到地址哈。），而且这个时候要明确你的参考基因组版本了！！！ b36/b37/hg18/hg19/hg38，记住b37和hg19并不是完全一样的，有些微区别哦！！！
Continue reading →

七 03

毕业生入深户完全指南

Posted on by

第一步：网上个人测评

申请人登录深圳市人力资源保障局官方网站(www.szhrss.gov.cn)，进入“网上办事”--“网上申办”--“深圳市人才引进(毕业生、在职人才引进)测评与申报系统”，注册个人账户，注册成功后通过个人用户登录系统选择 “毕业生接收”，根据系统提示填写个人信息，填报完成后，点击“保存”--点击“按当前填报信息测评”，系统将判断所填报人员是否符合毕业生接收政策并列出符合的政策条款。

也可以直接去测评网址，注册之后填一些信息https://sz12333.gov.cn/rcyj/

Ps:信息填写要真实，填写完了之后等待审核，一般三到五个工作日即可审核完毕，没什么特殊情况都会通过的，如果查看到自己审核通过了就可以进行第二步啦！

 

第二步：上门签订人事代理协议

符合毕业生接收政策的，即可与市人力资源局认可的人力资源代理机构签订个人申办委托办理协议，委托其办理毕业生接收手续。

上门需要带一些必备的资料，如下所示：

 

 

序号 材料名称
2 接收高等院校应届毕业生呈报表（收原件）
3 毕业生推荐表、成绩单（收原件）
4 学历及学位证书（申报时已毕业的验原件,收复印件；申报时未毕业的报到时验原件，不收复印件）
5 身份证（收复印件）
1
户口簿（户籍证明）

 

 

 

以上所有能带原件的都带上，然后所有原件都有复印一份！

代理机构有很多，大家选择自己最方便的， 我去的是深圳市人才交流服务中心（高新区分部） 

这个步骤需要上门，而且还需要排队，很可能需要排队两个到三个小时。还需要交钱，可能是260左右，可以刷卡。

PS：这个步骤因为要请假，所以大家一定要带全资料！！！办理很简单，主要是排队时间太长，办理完了会给你一个回执，你按照回执的提示15个工作日左右即可查看自己是否办理成功！如果成功了就再来一次，拿接收函！！！

 

第三步：用深圳市的接收函在学校拿报到证和户口迁移证

如果你是刚毕业，报到证还没有，那么这一步很简单，委托学校的同学帮忙即可。

如果你已经被开过报到证了（一般是遣返回老家啦），你就需要改派报到证啦！这个改派其实很简单，你需要自己看看你们学校改派流程，委托同学把新的报到证寄给你即可，如果你的档案还在学校就要求学校档案馆把你的档案通过机要传给深圳（15天左右），如果你的档案被遣返回家或者异地，那么你就要打电话去你的档案所在地要求他们帮你把你的档案通过机要传给深圳（15天左右）！

如果你的户口在学校，那么很简单，去你学校弄一个户口迁移证即可。

如果你的户口在老家，那么就麻烦了，还需要农转非什么的，看看你家里人的关系吧！

Ps：用深圳的接收函回学校成功拿到报到证和户口迁移证之后要随时上网查看自己的档案是否到达深圳。

 

第四步：拿介绍信和深圳市入户人员信息卡

这一个步骤不需要排队，在罗湖人才市场，需要身份证，毕业证，学位证，学历验证报告，报到证和户口迁移证原件及复印件各一份，缺一不可！！！

Ps：学历验证报告在学信网即可弄，请保证有效期至少一年以上！！！

第五步：去派出所办理户口身份证

这个需要预约！

这个需要预约！

这个需要预约！

重要的事情说三遍！如果你预约好了，那么你从罗湖人才市场拿到了介绍信和深圳市入户人员信息卡后就可以直接去派出所啦！！！但是如果你没有预约，你就得再等一个星期等到拿到预约时间后才能去派出所办理！

除了需要你在罗湖人才市场拿到了介绍信和深圳市入户人员信息卡，还需要数码照相回执和身份证，以及它们的复印件！！！

Ps：如果你是落户到高新园区派出所，那么你还有个近路，直接去迈瑞警务室也能完成落户流程！

到这里，落户就完成啦！十个工作日之后去派出所拿新的身份证即可！是不是非常简单呀小朋友们！

当然别忘了最后的彩蛋！网上深圳市新引进人才租房补贴系统

https://sz12333.gov.cn/szhr_pubtalent/talent_login.jsp

点击进入，有惊喜。

科未满30周岁、硕士未满35周岁、博士未满40周岁。租房补贴标准为：本科6000/人，硕士9000元/人；博士12000元/人。

 

总结一下：你需要请假三次或者四次，分别是去签人才引进代理协议，再去签人才引进代理协议的地方拿接收函，去罗湖人才市场拿介绍信和入户信息卡，去派出所办理落户及新身份证！

这个流程如果你仔细看了，而且保证按照流程走，当然，以你在各个单位拿到的最新资料为准，记住，各种材料宁可多带，也不能缺，一旦你少带了什么，没有人会跟你讲人情的，一切推倒重来！应该还算是蛮简单的，如果有任何疑问，欢迎咨询我QQ1227278128