#!/usr/bin/perl use strict; use warnings; my $sam=shift; my $vaf=$sam.".ps.tab"; my $roh=$sam."_ROHTableCall.bed"; my (%dat,$num); open my $rohfh, '<', $roh or die "Cannot open file: $roh $!\n"; while(<$rohfh>){ chomp; my ($chr, $start, $end, $flag, $snp)=split /\t/; $chr=~s/chr//; if($snp>150 or ($end-$start > 5000000 and $snp>50 ) ){ $num++; $dat{$chr}{$start.'_'.$end}='ROH'.$num; } } open my $vaffh, '<', $vaf or die "Cannot open file:$vaf $!\n"; open my $wfh, '>' , $vaf.".xls" or die "Cannot write file:$! \n"; $num=0; while(<$vaffh>){ chomp; next if /^#/; if(/^contig/){ print $wfh $_,"\tVAFadj\tROH\n"; next; } my ($chr,$pos, $nref, $nalt, $other, $af)=split /\t/; $chr=~s/chr//; next if $nalt+$nref<40; my $vaf=$nalt/($nalt+$nref); $vaf=$vaf >0.5 ? 1-$vaf : $vaf; $vaf=sprintf("%.3f", $vaf); my $roh; for my $reg (keys %{$dat{$chr} }){ my ($seg_start,$seg_end)=split /_/,$reg; if($pos>= $seg_start and $pos <= $seg_end){ $roh=$dat{$chr}{$reg}; ($num)=$roh=~/ROH(\d+)$/; last; } } if (defined $roh){ print $wfh $_,"\t", $vaf,"\t", $roh,"\n"; }else{ print $wfh $_,"\t", $vaf, "\t", 'HET'.$num,"\n"; } }