primer漏配问题解决

时间:2023-03-08 19:10:21

在对之前的ITS数据(454数据)做split时,发现有一些reads没有被匹配上,但是barcode能够完全匹配,虽然之后的primer在中间漏了一个碱基,导致后面的碱基全部误匹配,从而导致这条reads没有被匹配上的问题。

终于解决Qiime的问题后,使用 split_libraries.py 做切分,发现同样有这样的问题,Qiime并没有解决漏匹配的问题。

考虑如果用正常方法去做的话,对较小的异常数据需要花费N倍于正常数据的计算资源(包括硬件资源和运行时间),对于这个问题来说是非常不明智的。

由于primer长度在20左右,我的解决办法是,取末端6个碱基,在reads的这6个碱基相应位置左移1-2位,作为漏匹配1-2位的替代处理,这样既解决了漏匹配的问题,而且还能够使原正常的数据匹配速度进一步加快。

 #!/usr/bin/perl
use strict; my $usage = "usage:\nsplit.pl\tmapfile\t.fa_file\toutprefix\n";
die $usage unless @ARGV==; my $mapfile = shift @ARGV;
my $fafile = shift @ARGV;
my $outprefix = shift @ARGV; my %barcode;my $barcode_length = ;
my %primer;my $primer_length = ;
open MAP,$mapfile or die $!;
while(<MAP>){
chomp;
next if /^#/;
my @a = split /\s+/;
$barcode{$a[]} = $a[];
$barcode_length = length($a[]) unless $barcode_length;
die "barcode length do not match!" unless ($barcode_length == length($a[]));
$primer{$a[]} = $a[];
$primer_length = length($a[]) unless $primer_length;
die "primer length do not match!" unless ($primer_length == length($a[])); print "$barcode_length\t$primer_length\n";
}
close MAP; my %fa;
open FA,$fafile or die $!;
$/ = ">";
<FA>;
while(<FA>){
chomp;
my @a = split /\n/;
my $id = shift @a;
my $seq = join ("",@a);
@a = split (/\s+/,$id);
$id = shift @a;
$fa{$id} = $seq;
}
$/ = "\n";
close FA; open OUT,">$outprefix.fna" or die $!;
foreach my $id (sort keys %fa){
foreach my $sample (sort keys %barcode){
my $seq = substr($fa{$id},,$barcode_length);# print "$seq\n";
if ($barcode{$sample} eq $seq){
# print "barcode matched\n";
my $pri0 = substr($fa{$id},$barcode_length+$primer_length-,);
my $pri1 = substr($fa{$id},$barcode_length+$primer_length-,);
my $pri2 = substr($fa{$id},$barcode_length+$primer_length-,);
my $pri = substr($primer{$sample},$primer_length-,);
if ($pri0 eq $pri){
my $s = substr($fa{$id},$barcode_length+$primer_length,length($fa{$id})-$barcode_length-$primer_length);
print OUT ">$sample\t$id\n$s\n";
last;
}
if($pri1 eq $pri){ my $s = substr($fa{$id},$barcode_length+$primer_length-,length($fa{$id})-$barcode_length-$primer_length+);
print OUT ">$sample\t$id\n$s\n";
last;
}
if($pri2 eq $pri){
my $s = substr($fa{$id},$barcode_length+$primer_length-,length($fa{$id})-$barcode_length-$primer_length+);
print OUT ">$sample\t$id\n$s\n";
last;
}
}
}
}
close OUT;