tools/genotyping/pl/tmapfmt.pl

   1 #!/bin/env perl
   2 use strict;
   3 use warnings;
   4 use lib '/nas/RD_09C/resequencing/soft/lib';
   5 #use GalaxyXS::ChromByte;
   6 use Data::Dump qw(ddx);
   7
   8 unless (@ARGV > 0) {
   9     print "perl $0 <In_file> <Out_prefix> <ChrID>\n";
  10     exit 0;
  11 }
  12 =pod
  13 my %varToGT=(
  14         0   => 'A',
  15         1   => 'B',
  16         '-' => '-',
  17         2   => 'H',
  18 );
  19 =cut
  20 my ($Count,$population)=(0);
  21 my (%DAT,%deDup);
  22 open I,'<',$ARGV[0] or die $!;
  23 my $t=<I>;
  24 $t=~s/^#//;
  25 while (<I>) {
  26     chomp;
  27     my ($Pos,$D)=split /\s+/,$_,2;
  28         $DAT{$Pos}=$D;
  29         push @{$deDup{$D}},$Pos;
  30         #++$Count;
  31 }
  32 close I;
  33 no warnings;
  34 $population=split(/\s+/,(values %DAT)[0]);
  35 use warnings;
  36
  37 #ddx %deDup;
  38 #warn $population;
  39 #die;
  40
  41 my %TtoMarkers;
  42
  43 for my $D (keys %deDup) {
  44         my @T;
  45         @T=sort {$a<=>$b} @{$deDup{$D}};
  46         if (@T>1) {
  47                 $deDup{$D}=[ $T[0],$T[1] ];     # Remember the old version is $T[0,1] ! The correct is 0,-1 .
  48         } else {
  49                 $deDup{$D}=[ $T[0] ];
  50         }
  51         #$deDup{$D}=$max?[$min,$max]:[$min];
  52         $TtoMarkers{$D}=\@T;
  53 }
  54 my @OrderD=sort {$deDup{$a}->[0] <=> $deDup{$b}->[0]} keys %deDup;
  55 #ddx %deDup;
  56 my ($lastPos,$CurrPos)=0;
  57 for my $D (@OrderD) {
  58         $CurrPos = ($#{$deDup{$D}}==1)?($deDup{$D}->[1]):($deDup{$D}->[0]);
  59         warn "[!]$lastPos covers [",join(',',@{$deDup{$D}}),"].\n" if $CurrPos < $lastPos;
  60         $lastPos=$CurrPos;
  61 }
  62 $Count=keys %deDup;
  63
  64 open O,'>',$ARGV[1].'.intercross' or die $!;
  65 open OO,'>',$ARGV[1].'.order' or die $!;
  66 open G,'>',$ARGV[1].'.group' or die $!;
  67 print O "data type f2 intercross\n$population $Count\n";
  68 #for my $pos (sort {$a<=>$b} keys %DAT) {
  69 for my $D (@OrderD) {
  70         my $pos=join('-',@{$deDup{$D}});
  71         my $str=$D;
  72         $str =~ tr/012/ABH/;
  73         print O $ARGV[2],"m$pos\t",join("\t",split /\s+/,$str),"\n";
  74         print OO '*',$ARGV[2],"m$pos\n";
  75         print G substr($ARGV[2].'m'.$pos,0,15),'=',$ARGV[2],"m$pos\n";
  76         print G $ARGV[2],'m',$_,"\n" for @{$TtoMarkers{$D}};
  77         print G ">\n";
  78 }
  79 close G;
  80 close OO;
  81 close O;
  82
  83 __END__
  84 ./tmapfmt.pl ./dat2/Chr01_num Chr01.ic Chr01
  85 cat chrorder | while read a; do ./tmapfmt.pl ./dat2/${a}_num de$a $a;done &
  86
  87 ./tmap-1.0/pedmerge ril1206.Chr01.intercross tt
  88 cat chrorder | while read a; do ./tmap-1.0/pedmerge ./marker/de$a.intercross ./marker/de$a.pedigree;done &
  89 (Name beginning with Chr12m21531432- too long, truncating)
  90
  91 $ cat doim2.sh
  92 #!/bin/sh
  93 #$ -v PERL5LIB,PATH,PYTHONPATH,LD_LIBRARY_PATH,BLASTDB
  94 #$ -cwd -r y -l vf=3g
  95 #$ -S /bin/bash -t 1-12
  96 #$ -o /dev/null -e /dev/null
  97 SEEDFILE=./chrorder
  98 a=$(sed -n -e "$SGE_TASK_ID p" $SEEDFILE)
  99 eval ./tmap-1.0/tmap -i ./marker/de$a.pedigree ./marker/de$a.order > ./marker/de$a.imo 2> ./marker/de$a.imo.log