如果我有坐标,如何使用UCSC的Perl脚本提取DNA序列？

Question

如果我有坐标,如何使用基因组浏览器(UCSC)的Perl脚本提取DNA序列？

Answer 1

您可以将DAS序列请求传递给Perl脚本,该脚本解析包含序列的XML元素.

例如,以下是curlUCSC的DAS服务器的请求,丢弃标准错误,管道传递给parseSeq.pl:

$ curl http://genome.ucsc.edu/cgi-bin/das/hg19/dna?segment=1:10000,10999 2>/dev/null | parseSeq.pl

输出curl将是包含来自hg19人类基因组装配的1000碱基DNA序列的XML文档.请求从第一个染色体请求基数10000到10999(记住UCSC 从0开始).XML将包含一些对日志记录和错误检查有用的其他内容.

将XML 传递到Perl脚本之后,您可以使用Perl的XML :: Simple模块快速解析您想要的东西.

为了帮助您入门,您的parseSeq.pl文件可能从以下开始:

#!/usr/bin/perl -w                                                                                                                                                                                                                          

use strict;                                                                                                                                                                                                                                 
use XML::Simple;                                                                                                                                                                                                                            
use Data::Dumper;                                                                                                                                                                                                                           

my $xml = new XML::Simple;                                                                                                                                                                                                                  
my $ref = $xml->XMLin('-');                                                                                                                                                                                                                       

print Dumper $ref;

这个输出应该给你足够的开始拉动DNA序列$ref.