Hi
I"m still really new to BASH (and informatics in fact). I have extracted a contig into a .txt file, but I also have the information as one fasta within a multiple fasta file. I want to get one area of sequence that is in one contig "out" so I can BLAST it and do other things. What's a good way of doing this using BASH? (or something else simple in the command line).
So I want to...
take a sequence with about 750 000 bp in it, and I want only the sequence from position 60 000 to 75 000.
I feel like I should be able to combine grep and then control the grep command for specific lines into the file, then print to a .txt file...?
Thanks
Loren
1) If the sequence is available on NCBI, why not retrieve the fragment you want.
2) Blast the bigger fragment against any reference sequence (smaller) and download the aligned fragment.
3) grep reads the header, so what you are thinking might not be possible
4) trimming of sequence in sequence editor could be done (Aliview)
5) If you have only one or few sequences, text editor might also be an option (Notepad++ or Atom)
6) why not annotate and use the annotations ?
Hi, Loren
You can try bedtools getfasta.
Fistly, write the regeions (60 000 to 75 000) you want to extract into a BED/GFF/VCF file. Then use bedtools getfasta to extract the sequence you want.
This is what Bioperl is for:
perl -MBio::SeqIO -e 'my $sio=new Bio::SeqIO(-file=shift, -format=>"fasta"); my $seq=$sio->next_seq; print $seq->subseq(60 000, 75 000)."\n";' file.fa
Have not tested the syntax, but I think it should work...
#sed command that prints lines 60k-75k from file fasta.txt to out.txt
sed -n 60000,75000p fasta.txt > out.txt
- Badges
- Science topic
- Similar topics
- Filing