Hi! I tried this using a tab separated file called samples.txt and I'm getting the following error: NameError: name 'file' is not defined
I'm not sure what's wrong?
Hi,
I have a collection of 200 bacterial genome sequences with each file corresponding to separate bacterial species. The names of the sequence files are the corresponding bacterial species (e.g Nostoc punctiforme.fasta and 199 more). The sequences are downloaded from multiple databases so they have different header formats.
I am working on some perl scripts which can append my identifier ( e.g.NOPU ) to each fasta entry of the corresponding genome (So that each entry of the Nostoc punctiforme.fasta has NOPU appended after >) i am able to append and modify long identifiers through these scripts which asks for the file name and the identifier . In order to automate this process i have created a text file with two column, the first column for Filename and second column for identifier and each line having new entry e.g
**Filename** **My identifier**
Nostoc punctiforme NOPU
Prochlorococcus marinus PRMA
and 198 more entries
the idea behind this is to write a script which asks for a file name, searches the file name as well as identifier from the text file and appends the corresponding identifier to each entry of the fasta file and makes a new output file with the same name as the identifier or a script that asks for a directory of fasta file and picks fasta files one by one, look up the text file for identifier and do the needful.
However i am failing to have any lead in this direction and request people out there to help me.
Thanks
This is a common bioinformatics problem which we can generalise as "I want to read a fasta file and do something with it." Since you are using Perl, you should be using Bioperl and in particular, the SeqIO module.
You need to break the problem down into a series of steps. They are:
I'll assume that your fasta files and the file of file names are in the same directory. I'll also assume that your fasta file names are as given, i.e. they contain a space and have no suffix. This is not good practice and I'd suggest renaming the files to something like Nostoc_punctiforme.fa. You don't say whether the fasta files contain protein or nucleotide sequences but a good convention for suffixes is that used by the NCBI: .faa for protein, .ffn for nucleotide and .fna for nucleotide single sequence (e.g. a chromosome).
I'll also assume that (1) the columns of your text file are tab-separated, (2) the file is called ids.txt, (3) you want to write to a new file with the suffix .mod, containing fasta sequences with the modified header and (4) the identifier should be inserted directly after the ">" in the header.
#!/usr/bin/perl -w
use strict;
use Bio::SeqIO;
# open text file
my $file = "ids.txt";
open IN, $file;
while(<IN>) {
# remove line endings and split on tab
chomp;
my @line = split("\t", $_);
# check fasta file exists and open
if(-e $line[0]) {
my $fa = Bio::SeqIO->new(-file => $line[0], -format => "fasta");
# create new fasta file for modified output
my $outseq = Bio::SeqIO->new(-file => ">$line[0].mod", -format => "fasta");
while(my $seq = $fa->next_seq) {
my $id = $seq->display_id;
# modify header
$id = $line[1].$id;
# modified sequence object
my $newseq = Bio::Seq->new(-id => $id, -seq => $seq->seq);
# write modified sequence to file
$outseq->write_seq($newseq);
}
}
}
close IN;
Result using a small test file named Nostoc_punctiforme:
INPUT Nostoc_punctiforme
>seq1
agagagagagagaga
>seq2
gatacatagatacag
OUTPUT Nostoc_punctiforme.mod
>NOPUseq1
agagagagagagaga
>NOPUseq2
gatacatagatacag
One the solution would be to:
In python would look like so:
for sample in file('samples.txt'):
fname, postfix = sample.split()
output = file( "%s.new" % fname, 'wt')
for line in file(fname):
line = line.strip()
if line.startswith(">"):
line = line + postfix
output.write("%s\n" % line)
You can also invert the situation to first preload your samples into a dictionary like so:
store = {}
for sample in file('samples.txt'):
fname, postfix = sample.split()
store[fname]=postfix
fname = "one.fasta"
postfix = store[fname]
output = file( "%s.new" % fname, 'wt')
for line in file(fname):
line = line.strip()
if line.startswith(">"):
line = line + postfix
output.write("%s\n" % line)
Hi! I tried this using a tab separated file called samples.txt and I'm getting the following error: NameError: name 'file' is not defined
I'm not sure what's wrong?
try this: with "Filename My identifier" file (asuming tab delimited)
cat Filenames.txt | while read LINE; do FNAME="$(echo $LINE | cut -d"\t" -f1 )"; ID="$(echo $LINE | cut -d"\t" -f2)"; awk -v code="ID" '{sub(/^>/, "$code", $0); print}' ${FNAME}.fasta ; done
I haven't tested though
Log in to answer this question.