BOL: Owner

Install BLAST in Ubuntu/Linux and Window !

Poonam Mahapatra — Mon, 15 Mar 2021 12:35:50 -0500

#On ubuntu
sudo apt-get install ncbi-blast+

#Ubuntu Conda installation 
conda install -c bioconda blast

#Windows installation
First:
Download ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/2.7.1/ncbi-blast-2.7.1+-win64.exe with your web browser. 

#If that FTP links does not work, try HTTP instead http://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/2.7.1/ncbi-blast-2.7.1+-win64.exe

#Run this 
installer using the defaults, it should put BLAST under C:\Program Files\NCBI\blast-2.7.1+

Perl6 script to count ATGC !

Poonam Mahapatra — Tue, 28 Jan 2020 05:45:48 -0600

use v6;

my $default-input = "AGCTTTTCATTCTGACTGCAACGGGCAATATGTCTCTGTGTGGATTAAAAAAAGAGTGTCTGATAGCAGC";

sub MAIN($input = $default-input) {
    "{.map({ +$input.comb(/$_/) })}".say;
}

#I love perl v6

Perl script to convert fastq to fasta file

Poonam Mahapatra — Wed, 07 Mar 2018 04:23:26 -0600

#!/usr/bin/env perl

use strict;
use warnings;
use Bio::Factory::EMBOSS;

my $usage   = "perl $0 in.fq out.fa";
my $infile  = shift or die $usage;
my $outfile = shift or die $usage;

my $factory = Bio::Factory::EMBOSS->new;
my $seqret  = $factory->program('seqret'); # $seqret is a Bio::Tools::Run::EMBOSSApplication object

$seqret->run({-sequence => $infile,
              -sformat1 => 'fastq',
              -outseq   => $outfile,
              -osformat => 'fasta'});

Perl script to remove fasta sequences in multifasta file with certain length threshold

Poonam Mahapatra — Wed, 07 Mar 2018 04:16:32 -0600

#!/usr/bin/perl
use strict;
use warnings;
 
my $minlen = shift or die "Error: `minlen` parameter not provided\n";
{
    local $/=">";
    while(<>) {
        chomp;
        next unless /\w/;
        s/>$//gs;
        my @chunk = split /\n/;
        my $header = shift @chunk;
        my $seqlen = length join "", @chunk;
        print ">$_" if($seqlen >= $minlen);
    }
    local $/="\n";
}

Count the number of N in fasta file with Perl

Poonam Mahapatra — Fri, 04 Aug 2017 09:22:14 -0500

#!/usr/bin/perl

my ($h, $n, $l);

open(I,$ARGV[0]) or die($!);
while(){
  chomp;
  if(/^>/){
    $h=substr($_,1);
  }else{
    $n=($_=~tr/nN/nN/);
    $l=length($_);
    print $h,"\t",$l,"\t",$n,"\t",$n/($l-$n),"\n";
  }
}
close(I);

Perl script to calculate N50

Poonam Mahapatra — Fri, 09 Dec 2016 04:14:04 -0600
#! /usr/local/bin/perl5.8.3 -w ### ### Script Name : calcN50.pl ### Description : This script calculates the N50 of a collection of sequences in a fasta file ### Input : A fasta sequence file and the minimum length of sequence to consider it for the N50 calculation ### Output : None. The N50 value as well as the number of sequences >= this value is sent to STDOUT ### use strict; use Getopt::Long; use Bio::SeqIO; my $fasta_file; my $filter_length; GetOptions("fasta_file=s" => \$fasta_file, "length_filter=i" => \$filter_length); die "Usage: $0 [options] \t\t-f /path/to/fasta/file [Required] \t\t-l filter seq less than length (0 if no filter is to be done)\n\n" if !$fasta_file || !defined($filter_length) || $filter_length !~ /\d+/; ### ### Get the length of each sequence: ### my $in = Bio::SeqIO->new(-file => $fasta_file, -format => 'fasta'); my @seq_lengths = (); while (my $seq_obj = $in->next_seq){ next if length($seq_obj->seq) < $filter_length; push @seq_lengths, length($seq_obj->seq); } ### ### Determine the N50: ### my @vals; my $sum = 0; for my $val (@seq_lengths){ chomp $val; $sum += $val; push(@vals, $val); } my @sorted = sort { $a <=> $b } @vals; my $runningSum = 0; my $N50 = 0; foreach my $v (@sorted) { $runningSum += $v; if($runningSum >= ($sum / 2)) { $N50 = $v; last; } } ### ### Determine the number of sequences >= than the N50 in length: ### my $count = 0; foreach my $v (@sorted) { if($v >= $N50) { $count++; } } printf("N50: $N50 Count >= N50: $count\n");