<commit_before>#!/usr/bin/perl

=head1 NAME

 introgres_stats.pl
  print introgression start, end, and number of gene models within that region

=cut

=head1 SYPNOSIS

 introgres_stats.pl [-h] -i <inputtabfile> -a <gff3 annotation file>

=head2 I<Flags:>

=over


=item -i

<input_tab_file>           tabular file (mandatory)

=item -a

B<gff3_annotation_file>      gff3 file with gene models for counting number of genes in each region


=item -o 

B<out_file>             Output file for writing. Defaults to infile.sum

=item -h

B<help>                   print the help

=back

=cut

=head1 AUTHORS

  <NAME> <<EMAIL>>

=cut



use strict;
use warnings;

use Getopt::Long;
use File::Slurp qw (write_file);
use Pod::Usage;
use List::MoreUtils qw( minmax );

my ( $in_file, $gff, $out, $help); 

GetOptions (
    "infile|i=s"  => \$in_file,
    "gff|g|a=s"   => \$gff,
    "out|o=s"     => \$out,
    "help"        => \$help)   # flag
    or  pod2usage(-verbose  => 2);

if ($help || !$in_file || !$gff  || !$out)  { pod2usage(-verbose  => 2); } 

# clear existing out file 
write_file( $out , '');

open (IN, "<$in_file") || die "Cannot open input file $in_file.\n";

my $err = $out . ".err" ;
my %starts;
my %ends;

#parse the gff3 file 
open ( GFF3 , "<$gff" ) || die "Cannot open gff3 file $gff \n";
my %gff_genes;

while ( <GFF3> ) { 
    my $line = $_;
    chomp($line) ;
    my @fields = split ("\t" , $line ) ;
    my $chr = $fields[0];
    my $type = $fields[2];
    no  warnings "uninitialized";
    if ( $type eq "gene" ) {
	push @{ $starts{ $chr } } , $fields[3];
	push @{ $ends{ $chr } } , $fields[4];
    }
}

my %stats;

while ( <IN> ) { 
    my $line = $_;
    chomp($line) ; 
    my @fields = split ("\t" , $line) ; 
    my $num = $fields[0];
    my $chr = $fields[1];
    my $pos = $fields[2];
    my $snps = $fields[3];
    $stats{$num}->{chr} = $chr;
    push ( @ { $stats{$num}->{pos} } , $pos ); # {pos} is a list of positions, later the start and end will be extracted
    $stats{$num}->{snps} += $snps;
}


foreach my $introg_num (sort { $a <=> $b } keys %stats ) {
    my @positions = @{$stats{$introg_num}->{pos} };
    my ($start, $end )  = minmax @positions;
    my $size = scalar( @positions ) ;
    my $real_start = $start - ( ($end - $start)/($size-1)) + 1 ; #subtract bin size from first position
    my $bin_size =  ( ( $end - $start ) / ( $size -1 ) ) * ( $size );

    my $chr = $stats{$introg_num}->{chr};
    my $snps = $stats{$introg_num}->{snps};
    my $gene_count;
    my @filtered_e =  ( grep { $_ < $end } @{ $starts{ $chr } } );
    my @filtered_s =  grep { $_ > $start } @filtered_e;
    $gene_count = scalar( @filtered_s) ;
    write_file( $out, { append => 1 } , join("\t" , ($introg_num, $chr, $bin_size, $real_start, $end, $snps, $gene_count, "\n") ) ) ;
}<commit_msg>option or printing a gff file with the introgression regions. Useful for loading into jbrowse<commit_after>#!/usr/bin/perl

=head1 NAME

 introgres_stats.pl
  print introgression start, end, and number of gene models within that region

=cut

=head1 SYPNOSIS

 introgres_stats.pl [-h] -i <inputtabfile> -a <gff3 annotation file>

=head2 I<Flags:>

=over


=item -i

<input_tab_file>           tabular file (mandatory)

=item -a

B<gff3_annotation_file>      gff3 file with gene models for counting number of genes in each region


=item -o 

B<out_file>             Output file for writing. Defaults to infile.sum

=item --gff

B<gff3 out file>     Output a gff3 file

=item -p | --prefix

B<sample prefix>     the name of the sample to be used as prefix for introgression name in the gff3 output file
 
=item -h

B<help>                   print the help

=back

=cut

=head1 AUTHORS

  <NAME> <<EMAIL>>

=cut



use strict;
use warnings;

use Getopt::Long;
use File::Slurp qw (write_file);
use Pod::Usage;
use List::MoreUtils qw( minmax );

my ( $in_file, $gff, $out, $gff_out, $prefix, $help); 

GetOptions (
    "infile|i=s"  => \$in_file,
    "annot|a=s"   => \$gff,
    "out|o=s"     => \$out,
    "gff|g=s"       => \$gff_out,
    "p|prefix=s"  => \$prefix,
    "help"        => \$help)   # flag
    or  pod2usage(-verbose  => 2);

if ($help || !$in_file || !$gff  || !$out)  { pod2usage(-verbose  => 2); } 

# clear existing out file 
write_file( $out , '');

open (IN, "<$in_file") || die "Cannot open input file $in_file.\n";

my $err = $out . ".err" ;
my %starts;
my %ends;

#parse the gff3 file 
open ( GFF3 , "<$gff" ) || die "Cannot open gff3 file $gff \n";
my %gff_genes;

while ( <GFF3> ) { 
    my $line = $_;
    chomp($line) ;
    my @fields = split ("\t" , $line ) ;
    my $chr = $fields[0];
    my $type = $fields[2];
    no  warnings "uninitialized";
    if ( $type eq "gene" ) {
	push @{ $starts{ $chr } } , $fields[3];
	push @{ $ends{ $chr } } , $fields[4];
    }
}

my %stats;

while ( <IN> ) { 
    my $line = $_;
    chomp($line) ; 
    my @fields = split ("\t" , $line) ; 
    my $num = $fields[0];
    my $chr = $fields[1];
    my $pos = $fields[2];
    my $snps = $fields[3];
    $stats{$num}->{chr} = $chr;
    push ( @ { $stats{$num}->{pos} } , $pos ); # {pos} is a list of positions, later the start and end will be extracted
    $stats{$num}->{snps} += $snps;
}


foreach my $introg_num (sort { $a <=> $b } keys %stats ) {
    my @positions = @{$stats{$introg_num}->{pos} };
    my ($start, $end )  = minmax @positions;
    my $size = scalar( @positions ) ;
    my $real_start = $start - ( ($end - $start)/($size-1)) + 1 ; #subtract bin size from first position
    my $bin_size =  ( ( $end - $start ) / ( $size -1 ) ) * ( $size );

    my $chr = $stats{$introg_num}->{chr};
    my $snps = $stats{$introg_num}->{snps};
    my $gene_count;
    my @filtered_e =  ( grep { $_ < $end } @{ $starts{ $chr } } );
    my @filtered_s =  grep { $_ > $start } @filtered_e;
    $gene_count = scalar( @filtered_s) ;
    write_file( $out, { append => 1 } , join("\t" , ($introg_num, $chr, $bin_size, $real_start, $end, $snps, $gene_count, "\n") ) ) ;
    if ( $gff) { 
	write_file( $gff_out , { append => 1 }, join("\t" , ( $chr, "Introgression" , "match" , $real_start , $end, "." , "+" , "." , "ID=" . $prefix . "_introgression" . $introg_num , "\n")));
    }
}
{'fraction_non_alphanumeric': 0.18881006127055788, 'fraction_numerical': 0.009351821992905514, 'mean_word_length': 2.786935286935287, 'pattern_counts': {'":': 0, '<': 37, '<?xml version=': 0, '>': 62, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '4310595', 'n_tokens_mistral': 2548, 'n_tokens_neox': 2457, 'n_words': 683}