annotate egglib/CalculateDiversityIndexes.pl @ 14:d15869b3731a draft

Deleted selected files
author dereeper
date Tue, 08 Jan 2019 08:46:57 -0500
parents 98c37a5d67f4
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
9
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
1 #!/usr/bin/perl
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
2
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
3 use strict;
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
4
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
5 use Getopt::Long;
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
6 use Bio::SeqIO;
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
7
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
8
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
9 my $usage = qq~Usage:$0 <args> [<opts>]
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
10 where <args> are:
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
11 -i, --input <FASTA input>
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
12 -o, --output <output filename>
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
13 -d, --directory <directory of egglib package>
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
14 ~;
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
15 $usage .= "\n";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
16
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
17 my ($infile,$outfile,$dir_exe);
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
18
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
19
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
20 GetOptions(
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
21 "input=s" => \$infile,
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
22 "output=s" => \$outfile,
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
23 "directory=s"=> \$dir_exe
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
24 );
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
25
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
26
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
27 die $usage
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
28 if ( !$infile || !$outfile || !$dir_exe);
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
29
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
30
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
31 my $EGGSTATS_EXE = "$dir_exe/egglib-2.1.5/bin/eggstats";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
32
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
33 my %gene_alignments;
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
34 my $in = Bio::SeqIO->new(-file => $infile , '-format' => 'Fasta');
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
35 while ( my $seq = $in->next_seq() )
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
36 {
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
37 my $id = $seq -> id();
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
38 my $sequence = $seq -> seq();
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
39 my ($gene,$ind,$num_allele) = split("_",$id);
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
40 $gene_alignments{$gene}.= ">$id\n$sequence\n";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
41 }
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
42
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
43 open(OUT,">$outfile");
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
44 foreach my $gene(keys(%gene_alignments))
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
45 {
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
46 open(F,">$gene.egglib_input.fa");
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
47 print F $gene_alignments{$gene};
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
48 close(F);
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
49
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
50 my $results_egglib = `$EGGSTATS_EXE $gene.egglib_input.fa`;
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
51
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
52 # parse Seqlib output
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
53 if ($results_egglib)
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
54 {
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
55 my %egglig_stats;
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
56 my @eggstats = split(/^/,$results_egglib);
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
57 foreach my $eggstat(@eggstats)
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
58 {
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
59 my ($desc,$value) = split(/: /,$eggstat);
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
60 chomp($value);
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
61 $egglig_stats{$desc} = $value;
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
62 }
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
63 print OUT "$gene;";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
64 print OUT $egglig_stats{"Total number of sequences"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
65 print OUT $egglig_stats{"Total number of sites"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
66 print OUT $egglig_stats{"Number of analyzed sites"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
67 print OUT $egglig_stats{"S"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
68 print OUT $egglig_stats{"thetaW"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
69 print OUT $egglig_stats{"Pi"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
70 print OUT $egglig_stats{"D"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
71 print OUT $egglig_stats{"number of haplotypes"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
72 print OUT $egglig_stats{"haplotypes diversity"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
73 print OUT $egglig_stats{"Fay and Wu H"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
74 print OUT $egglig_stats{"Fst"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
75 print OUT $egglig_stats{"Snn"} . ";";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
76 print OUT "\n";
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
77 unlink("$gene.egglib_input.fa");
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
78 }
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
79 }
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
80 close(OUT);
98c37a5d67f4 Uploaded
dereeper
parents:
diff changeset
81