Hi all,

i m trying to export the results of a simple BLASTn query (16S ribosomal RNA) against Genbanks nr database using library 'annotate' and the function 'blastsequences' (see: http://www.bioconductor.org/packages/release/bioc/manuals/annotate/man/annotate.pdf)

Function works well, hitlist is retrieved correctly as dataframe or multiple alignment, but there is no way to write the results into Fasta format!(?) I found hints that this might be possible via 'Biostrings' (which is essential part of the function). Can anyone provide a solution for this? Help is much appreciated!

Best regards,

Benjamin

See example below:

library(Biostrings)
library(annotate)
blastSequences(x = "AGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACACATGCAAGTCGAACGGTAACAGGAAG
CAGCTTGCTGCTTTGCTGACGAGTGGCGGACGGGTGAGTAATGTCTGGGAAACTGCCTGATGGAGGGGGA
TAACTACTGGAAACGGTAGCTAATACCGCATAACGTCGCAAGCACAAAGAGGGGGACCTTAGGGCCTCTT
GCCATCGGATGTGCCCAGATGGGATTAGCTAGTAGGTGGGGTAACGGCTCACCTAGGCGACGATCCCTAG
CTGGTCTGAGAGGATGACCAGCAACACTGGAACTGAGACACGGTCCAGACTCCTACGGGAGGCAGCAGTG
GGGAATATTGCACAATGGGCGCAAGCCTGATGCAGCCATGCNGCGTGTATGAAGAAGGCCTTCGGGTTGT
AAAGTACTTTCAGCGGGGAGGAAGGGAGTAAAGTTAATACCTTTGCTCATTGACGTTACCCGCAGAAGAA
GCACCGGCTAACTCCGTGCCAGCAGCCGCGGTAATACGGAGGGTGCAAGCGTTAATCGGAATTACTGGGC
GTAAAGCGCACGCAGGCGGTTTGTTAAGTCAGATGTGAAATCCCCGGGCTCAACCTGGGAACTGCATCTG
ATACTGGCAAGCTTGAGTCTCGTAGAGGGGGGTAGAATTCCAGGTGTAGCGGTGAAATGCGTAGAGATCT
GGAGGAATACCGGTGGCGAAGGCGGCCCCCTGGACGAAGACTGACGCTCAGGTGCGAAAGCGTGGGGAGC
AAACAGGATTAGATACCCTGGTAGTCCACGCCGTAAACGATGTCGACTTGGAGGTTGTGCCCTTGAGGCG
TGGCTTCCGGANNTAACGCGTTAAGTCGACCGCCTGGGGAGTACGGCCGCAAGGTTAAAACTCAAATGAA
TTGACGGGGGCCGCACAAGCGGTGGAGCATGTGGTTTAATTCGATGCAACGCGAAGAACCTTACCTGGTC
TTGACATCCACGGAAGTTTTCAGAGATGAGAATGTGCCTTCGGGAACCGTGAGACAGGTGCTGCATGGCT
GTCGTCAGCTCGTGTTGTGAAATGTTGGGTTAAGTCCCGCAACGAGCGCAACCCTTATCCTTTGTTGCCA
GCGGTCCGGCCGGGAACTCAAAGGAGACTGCCAGTGATAAACTGGAGGAAGGTGGGGATGACGTCAAGTC
ATCATGGCCCTTACGACCAGGGCTACACACGTGCTACAATGGCGCATACAAAGAGAAGCGACCTCGCGAG
AGCAAGCGGACCTCATAAAGTGCGTCGTAGTCCGGATTGGAGTCTGCAACTCGACTCCATGAAGTCGGAA
TCGCTAGTAATCGTGGATCAGAATGCCACGGTGAATACGTTCCCGGGCCTTGTACACACCGCCCGTCACA
CCATGGGAGTGGGTTGCAAAAGAAGTAGGTAGCTTAACTTCGGGAGGGCG
" ,database = "nr", hitListSize="20", program = "blastn", as="data.frame")

 Hi Benjamin

 when you search for function it is good to list functions in the packages you already use  or type "??fasta "  or  search R.docs website.  Biostrings developers made already writeXStringSet
(In your example i used only first 70 NN as i hoped response from website would be faster ... )

library(annotate);
library(Biostrings) ;ls("package:Biostrings")

BLrez =blastSequences(x = "AGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACACATGCAAGTCGAACGGTAACAGGAAG" ,database = "nr", hitListSize="20", program = "blastn", as="data.frame")
 
BLfasta = DNAStringSet (BLrez$Hsp_hseq)
names (BLfasta) = as.character(BLrez$Hit_id)
 
writeXStringSet (BLfasta, file="rRNA.fasta.txt", format="fasta", width=80)

btw. great  R & sequencing tutorial (bit old)  is here .
Branko

Dear Branko,

thanks so much for your help, very much appreciated, that worked perfectly! Actually i already came across 'writeXStringSet' but struggled with its use. Thanks again!

Best regards,

Benjamin