Bioinformatics Programming
Extract FASTA sequences based on sequence length using Perl
Here are simple Perl scripts to filter out FASTA sequences from a multi-fasta file based on sequence length.
Let’s say our input file consisting of multiple FASTA sequences is ‘input.fasta’.
#!/usr/bin/perl
use strict;
use warnings;
my ($infile, $minlen) = @ARGV;
{
local $/=">";
while(<$infile>) {
chomp;
next unless /\w/;
my @keep = split /\n/;
my $header = shift @keep;
my $seqlen = length join "", @keep;
if($seqlen >= $minlen){
print ">$_";
}
}
local $/="\n";
}
exit;
Save this Perl script as ‘extractfasta.pl‘ and run in the terminal as
$ perl extractfasta.pl input.fasta <minlen> > output.fasta
For example,
$ perl extractfasta.pl input.fasta 100 > output.fasta
If you want to set a maximum length limit as well, then use the following script.
#!/usr/bin/perl
use strict;
use warnings;
my ($infile, $minlen, $maxlen) = @ARGV;
{
local $/=">";
while(<$infile>) {
chomp;
next unless /\w/;
my @keep = split /\n/;
my $header = shift @keep;
my $seqlen = length join "", @keep;
if($seqlen >= $minlen){
print ">$_";
}
}
local $/="\n";
}
exit;
Save this Perl script as ‘extractfasta.pl‘ and run in the terminal as
$ perl extractfasta.pl input.fasta <minlen> <maxlen> > output.fasta
For example,
$ perl extractfasta.pl input.fasta 100 350 > output.fasta
Fetching SMILE strings for a given number of SDF files of chemical compounds is not such a trivial task. We can quickly obtain them using RDKit or OpenBabel. But what if you don’t have SDF files of ligands in the first place? All you have is Ligand IDs from PDB. If they are a few then you can think of downloading SDF files manually but still, it seems time-consuming, especially when you have multiple compounds to work with. Therefore, we provide a Python script that will read all Ligand IDs and fetch their SDF files, and will finally convert them into SMILE strings. (more…)
Bioinformatics Programming
How to get secondary structure of multiple PDB files using DSSP in Python?
In this article, we will obtain the secondary structure of multiple PDB files present in a directory using DSSP [1]. You need to have DSSP installed on your system. (more…)
Bioinformatics Programming
vs_analysis_compound.py: Python script to search for binding affinities based on compound names.
Previously, we have provided the vs_analysis.py script to analyze virtual screening (VS) results obtained from Autodock Vina. In this article, we have provided another script to search for binding affinity associated with a compound. (more…)
[email protected]
April 13, 2021 at 3:52 am
I got an error when usong the second command:
readline() on unopened filehandle at extractfasta.pl line 7
this is what I run:
perl extractfasta.pl /path-to/BAC4A_L00M_R1_001.fasta 50 100 > 100_maxln.fasta
the input fasta is ok , dont know what is wrong ):