使用三代测序数据能获得较好的、甚至完整的基因组序列。通过检测基因组序列两端的属于端粒（Telomere）的特定的重复序列，可以知道基因组组装是否得到完整的染色体水平的序列。若在序列中间检测到端粒序列，可以知道基因组组装过程中对Contigs有错误的连接。
在染色体序列首尾存在端粒序列。在人类中，端粒序列由重复单元TTAGGG，串联重复约2500次组成。 不同的物种，端粒的重复单元可能不一样，可以在端粒数据库中查询。 我对一种子囊菌使用PacBio测序数据进行了基因组组装，得到12条序列，发现大部分序列首尾均出现了TTAGGG/CCCTAA的串联重复序列。我对一种昆虫物种PacBio组装基因组序列进行分析，发现序列首尾出现一些长度较长微卫星重复序列，而没有固定的重复单元，这和端粒数据库中的结果一致。
编写Perl程序对全基因组序列的端粒序列进行搜索，查看基因组序列的完整情况：
#!/usr/bin/perl
use strict;
use Getopt::Long;
my $usage = <<USAGE;
Usage:
$0 genome.fasta > telomere_info.txt
大部分物种端粒序列的重复单元是TTAGGG/CCCTAA。本程序能在基因组中寻找端粒重复单元的串联重复序列，并给出位点信息。
--split-length <int> default: 100000
--overlap-length <int> default: 10000
程序会将每条序列打断后进行重复单元搜索。这两个参数设置打断的序列长度和相邻两序列之间的重叠长度。
--repeat-unit <string> default: TTAGGG
设置重复单元碱基序列，该重复单元的反向互补也将作为重复单元进行搜索。可以在端粒数据库（http://telomerase.asu.edu/sequences_telomere.html）中寻找目标端粒重复单元。
vertebrate sp. TTAGGG
plants sp. TTTAGGG
Pezizomycotina TTAGGG
--min-repeat-num <int> default: 4
设置重复单元最小重复次数。
USAGE
if (@ARGV==0){die $usage}
my ($splitLength, $overlapLength, $repeatunit, $minRepeatNum);
GetOptions(
"split-length:i" => \$splitLength,
"overlap-length:i" => \$overlapLength,
"repeat-unit:s" => \$repeatunit,
"min-repeat-num:s" => \$minRepeatNum,
);
$splitLength ||= 100000;
$overlapLength ||= 10000;
$repeatunit ||= "TTAGGG";
$repeatunit = uc($repeatunit);
my $repeatunit_reverse = reverse $repeatunit;
$repeatunit_reverse =~ tr/ATCG/TAGC/;
$minRepeatNum ||= 4;
# 读取基因组序列
open IN, $ARGV[0] or die "Can not open file $ARGV[0], $!\n";
my (%seq, $seq_id);
while (<IN>) {
chomp;
if (m/^>(\S+)/) {
$seq_id = $1;
}
else {
$_ = uc($_);
$seq{$seq_id} .= $_;
}
}
close IN;
# 将基因组序列打断
my (%seq_split, %seq_length);
foreach my $id (keys %seq) {
my $seq = $seq{$id};
my $length = length($seq);
$seq_length{$id} = $length;
my $pos = 0;
while ($pos < $length) {
$seq_split{$id}{$pos} = substr($seq, $pos, $splitLength + $overlapLength);
$pos += $splitLength;
}
}
print "SeqID\tSeqLength\tStart\tEnd\tLength\tType\n";
foreach my $id (sort keys %seq_split) {
foreach my $pos (sort {$a <=> $b} keys %{$seq_split{$id}}) {
my $seq = $seq_split{$id}{$pos};
while ($seq =~ m/(($repeatunit){$minRepeatNum,})/g) {
my $length = length($1);
my $end = pos($seq);
$end = $end + $pos;
my $start = $end - $length + 1;
print "$id\t$seq_length{$id}\t$start\t$end\t$length\t$repeatunit\n";
}
while ($seq =~ m/(($repeatunit_reverse){$minRepeatNum,})/g) {
my $length = length($1);
my $end = pos($seq);
$end = $end + $pos;
my $start = $end - $length + 1;
print "$id\t$seq_length{$id}\t$start\t$end\t$length\t$repeatunit_reverse\n";
}
}
}