Thankyou so much..your script works perfectly. ;) Even though i've merged the multifasta sequence as single sequence,then my program also worked nicely.
i want to replace name of the sequences in fasta file after '>' sign with the protein id within the header line.i wrote a perl script but it only prints the firstline of the sequence itself(given below).
>XP_020088267.1
ATGGCTGATGCTGAGGATATTCAGCCACTCGTCTGTGACAATGGAACTGGAATGGTGAAGGCTGGATTTGCTGGTGATGA
>XP_020087759.1
ATGGTTGTTTCTAGCTCACCTAAAACTGCATATGATGCTTGGAGAATTATTGAGGCATATTTTCTTGATAAGACTGCTTC
>XP_020089204.1
ATGACGGCGACGACTTCTGACGCTATGGAGGAGGATCCGGCGCCGTCTTCTGATATGACGGAGGAGGAGGGGGGCGACCC
Here is the script;
my ($data,$seqline);
$infile=$ARGV[0];
open(IN,"$infile") || die "Can't open input file";
while($data=<IN>)
{
chomp $data;
$seqline=<IN>;
if($data=~s/^>//g)
{
my ($protn_id)= (split "=",$data)[4];
print ">$protn_id\n","$seqline\n";
}
}
close IN;
3 answers
I think the Error is in structures of input fasta file. Your script seems not work fine unless your fasta file is as follows, i guess.
>seq1 something=something=something=something=XP_020088267.1
ATGGCTGATGCTGAGGATATTCAGCCACTCGTCTGTGACAATGGAACTGGAATGGTGAAGGCTGGATTTGCTGGTGATGA
>seq2 something=something=something=something=XP_020087759.1
ATGGTTGTTTCTAGCTCACCTAAAACTGCATATGATGCTTGGAGAATTATTGAGGCATATTTTCTTGATAAGACTGCTTC
>seq3 something=something=something=something=XP_020089204.1
ATGACGGCGACGACTTCTGACGCTATGGAGGAGGATCCGGCGCCGTCTTCTGATATGACGGAGGAGGAGGGGGGCGACCC
Then the following script may be more robust.
my $infile=$ARGV[0];
open(IN,$infile) || die "Can't open input file";
while(my $line=<IN>){
$line =~ s/[\r\n]//g;
if($line=~ /^>/){
my ($protn_id)= (split "=",$line)[4];
print ">$protn_id\n";
}else{
print $line."\n";
}
}
close IN;
BTW, Perl v5 is an old language. Perl v6 is a very minor language.
If you are a beginner, I recommend python 3 as jrj.healey says.
That can be done using a Perl-one-liner:
$ perl -lne '(m/>.+\[protein_id=(.+?)\]/) ? print ">$1" : print $_' < seqs.fa
>XP_020098752.1
ATGGTTGCCACTAAGTTGATGATGACCTCTTTAATCTTAGTTCAACTGTGGGTGCTTATGCCACTGATGGCGTGTGGTAC
AACGTTAGATCCCATGAGAGAGAGGTATGAACAATGGATTAGCCGATATAGCCGAGTCTACAAGGATAAGAACGAGAAGG
AGTGGCGCTTTAGGATATATGAATCCAACGTCCAGCTCATCAACATCTTTAATACCATTAGTGAGGAGTACAAGCTCATT
GACAACAAGTTTGCCGACCTAACTAGTGAGGAGTTCAAGGCCAAGTCTGTTTGCTTAAGGGATCTCCGTAATCATCGTCC
GCCTTCTCGACAGTCGCAGCAGTAGAAGGCATCAACAAAATTAAGGCGGGTAGATTGGTAG
GTCAGTAGCAATAGATGCTGGGGGTTTCGCCTTCCAGTTCTACTCAAAGGGCATCTTCACC
GATTGCCATGAAGCCTTCCTATCCTCTCAAGATAGATTAA
>XP_020098740.1
ATGGCGGTGCGGGCAGCTGGTCTGCTGGCAGCGTTGGTTGTGGGTTTGGCCGCGGTGTACTGCGCGATGGACCCCCTACG
GCTGAGCGCCGTAGCCGACTTCCCGGGCTTCGAGAGCCATCCCGTAGAGCTTCCCCCTTGGTCGGAGCTGCCGGCCGCCA
GGGACGCCGAGGATCGGCTGCGGAGAGCGGAGATCCGCTTCCTGAACCAGGTGCAGGGCCCCGAGAGCATCGCCTTCGAC
GTACGGCCCAGAGGGGGAGTTGCTTGAAATTCTGGAAGACCGGCAGGGGAAGGTTGTTAGGGCAGTTAGCGAAGTCGAAG
AGAAGGATGGGAAGCTTTGGATAGGATCAGTGCTCATGCCATTTATTGCCGTTTATTGA
>XP_020084954.1
ATGCGAGCTCGTGCAGAGTCGATAATGGGCGGCGATCAGGGGAAGACAGCAATGGCGCAGAGGAAGCTTCTTCTCCGTGG
TCCGACGGCACTCGCCCAACGGCAGCCGGATGCCGCCTCTGCCTCAAAGCCACTGGGCCGGCGGCGGATAGCGGAGATGG
CTGGGGAGACGGCAGCGGAGTGTGCGGCCATCTGGTGCTGCTGTCCCTGCGGCCTCCTCAACCTCGTCGTCCTCGCCCTC
GTCAAACTTCCCGCCGGGCTCGTCATCCGCGCTCTGCGCCGCCAAAGGCGAGTGATGAGGAAGCACAGCAGAGGAGGAGC
CGCTGCCCTGAGGCTGCTGGGAGGGCCCAAAGGCAAAGGCAGAGGCAGAGGCTCAACTGGCGAAGGCGAAGGCGAAGGCG
Explanation:
perl -lneactivates perl in execute mode and new line read( condition ) ? exec1 : exec2this is a if-then-else abbreviated, if condition, then do exec1, otherwise do exec2m/>.+\[protein_id=(.+?)\]/this is a regular expression evaluation, looks for a line with >, then followed by any chars (.+) and looks for a string [protein_id=, the text after that will be recorded in in a match variable ($1), it is marked as any char but expanded until a ] is seen (.+?)].- Finally, if the match exists, it will print a new line as >$1, otherwise will print the line as read
Thank you.. Glad to know about the perl one-liner.. :)
You have an if($data=~s/^>//g), what happens when a line does not start at >? You are just missing an else and a print.
Log in to answer this question.
Then why don't you show your code?
I added code markup to your post for increased readability. You can do this by selecting the text and clicking the 101010 button. When you compose or edit a post that button is in your toolbar, see image below:
Thanks a lot... Please help me on my program
Is there a reason you aren’t just using
Bio::Perl?Hii.. I'm a beginner in Perl scripting and scripting in general.Also i'm not familiar with bioperl modules. that's why..
It's generally easier to use something that already exists, rather than reinventing the wheel.
I'm biased, but Perl is also not the easiest language to start with as a beginner (consider python).
Can you show us your input data too (just 3 or 4 sequences would do). It's not clear where the headers are coming from (I'm guessing they're already in the fasta and you're just deleting part of the header judging by your perl substitution.
Here with the example of my input file..
and the output should be like the following.. thanks.
An example of your input file and how do you want the output could be helpful.
BTW I proficient in Perl and Python, still using primarily Perl ;)