diff Admixture.pl @ 5:97c9c8daa3c3 draft

planemo upload
author dereeper
date Wed, 13 Apr 2016 07:51:13 -0400
parents 58df6910f1c3
children 3f66f32dc5d9
line wrap: on
line diff
--- a/Admixture.pl	Tue Apr 12 09:31:56 2016 -0400
+++ b/Admixture.pl	Wed Apr 13 07:51:13 2016 -0400
@@ -1,7 +1,9 @@
 #!/usr/bin/perl
 
 use strict;
+use Switch;
 use Getopt::Long;
+use Bio::SeqIO;
 use File::Basename;
 
 my $usage = qq~Usage:$0 <args> [<opts>]
@@ -11,10 +13,11 @@
     -k, --kmin          <K min. int>
     -m, --maxK          <K max. int>
     -d, --directory     <temporary directory>
+    -t, --threshold     <threshold admixture proportion for group assignation>
 ~;
 $usage .= "\n";
 
-my ($input,$output,$kmin,$kmax,$directory);
+my ($input,$output,$kmin,$kmax,$directory,$threshold);
 
 
 GetOptions(
@@ -23,6 +26,7 @@
 	"kmin=s"       => \$kmin,
 	"maxK=s"       => \$kmax,
 	"directory=s"  => \$directory,
+	"threshold"    => \$threshold,
 );
 
 
@@ -63,7 +67,7 @@
 my %errors;
 for (my $k = $kmin; $k <= $kmax; $k++)
 {
-	system("admixture --cv $input.bed $k >>$directory/log.$k 2>&1");
+	system("/apps/www/sniplay.cirad.fr/tools/admixture/admixture_linux-1.23/admixture --cv $input.bed $k >>$directory/log.$k 2>&1");
 	my $cv_error_line = `grep -h CV $directory/log.$k`;
 	if ($cv_error_line =~/: (\d+\.*\d*)$/)
 	{
@@ -73,6 +77,7 @@
 	system("echo '\n\n====================================\n\n' >>$directory/logs");
 
 	open(my $O2,">$basename.$k.final.Q");
+	open(my $O3,">$directory/groups.$k");
 	open(my $O,"$basename.$k.Q");
 	my %hash_groupes;
 	my %hash_indv;
@@ -88,7 +93,7 @@
 		my $ind = $individus[$i];
 		for (my $j = 0; $j <$k; $j++){
 			my $val = $infos[$j];
-			if ($val > 0.5){$group = "Q$j";}
+			if ($val > ($threshold/100)){$group = "Q$j";}
 		}
 		if ($ind){      
 			$hash_indv{$ind} = join("	",@infos);
@@ -102,12 +107,12 @@
 		my @inds = split(",",$hash_groupes{$group}{"ind"});
 		foreach my $ind(@inds){
 			if ($ind =~/\w+/){
-				#print $O3 "$ind;$group\n";
+				print $O3 "$ind;$group\n";
 				print $O2 $ind."	".$hash_indv{$ind}. "\n";
 			}
 		}
 	}
-	#close($O3);
+	close($O3);
 	close($O2);
 
 	system("cat $basename.$k.final.Q >>$directory/outputs.Q");
@@ -141,7 +146,7 @@
 close(BEST2);
 
 system("cp -rf $directory/log.$best_K $directory/log");
+system("cp -rf $directory/groups.$best_K $directory/groups");
 
 
 
-