annotate MDSplot/MDSbasedOnIBSmatrix.pl.org @ 1:420b57c3c185 draft

Uploaded
author dereeper
date Fri, 10 Jul 2015 04:39:30 -0400
parents
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
1
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
1 #!/usr/bin/perl
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
2
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
3 use strict;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
4 use Switch;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
5 use Getopt::Long;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
6 use Bio::SeqIO;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
7
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
8 my $PLINK_EXE= "/apps/www/sniplay.cirad.fr/tools/plink/plink-1.07-x86_64/plink";
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
9
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
10 my $usage = qq~Usage:$0 <args> [<opts>]
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
11 where <args> are:
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
12 -i, --in <input>
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
13 -o, --out <output>
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
14 ~;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
15 $usage .= "\n";
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
16
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
17 my ($in,$out);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
18
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
19
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
20 GetOptions(
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
21 "in=s" => \$in,
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
22 "out=s" => \$out
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
23 );
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
24
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
25 die $usage
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
26 if ( !$in || !$out);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
27
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
28
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
29 my $plink_command = $PLINK_EXE . " --file $in --noweb --cluster --matrix --mds-plot 2 --out $out >>$in.plink.log 2>&1";
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
30 system($plink_command);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
31
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
32 my $awk_cmd = "awk \{\'print \$1\'\} $in.ped";
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
33 my $inds = `$awk_cmd`;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
34 my @individuals = split("\n",$inds);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
35
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
36
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
37 open(my $OUT,">$out.mds_plot.txt");
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
38 my $go = 0;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
39 open(my $O,"$out.mds");
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
40 while(<$O>)
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
41 {
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
42 if ($go)
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
43 {
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
44 my $line = $_;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
45 $line =~s/\n//g;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
46 $line =~s/\r//g;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
47 my @i = split(/\s+/,$line);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
48 my $ind = $i[1];
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
49 print $OUT "$ind ".$i[4]." ".$i[5]."\n";
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
50 }
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
51 if (/C1/){$go = 1;}
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
52 }
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
53 close($O);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
54 close($OUT);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
55
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
56
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
57 my $j = 0;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
58 open(my $IBS,">$out.ibs_matrix.txt");
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
59 print $IBS "Individuals " . join("\t",@individuals)."\n";
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
60 open(my $O2,"$out.mibs");
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
61 while(<$O2>)
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
62 {
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
63 my $line = $_;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
64 $line =~s/\n//g;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
65 $line =~s/\r//g;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
66 my @i = split(/\s+/,$line);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
67 print $IBS $individuals[$j]. " ". join("\t",@i)."\n";
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
68 $j++;
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
69 }
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
70 close($O2);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
71 close($IBS);
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
72
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
73
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
74
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
75
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
76
420b57c3c185 Uploaded
dereeper
parents:
diff changeset
77