I want to identify the position and amino acid change in a specific region of interest in my sequences, and store that information in a table.

Is it possible to do something like this in R, by maybe using the bioconductor package? I have managed to do a simple alignment of the sequences with AlignSeqs() from the DECIPHER package, but I am unable to pull out the differences in the sequences automatically. I start with FASTA files.

I want to end up with something like this:

Isolate ID    Reference_AA   Sample_AA   Pos
1             S              T           254
2             T              D           200
3             L              A           230

I have a reference sequence which is 74 AA long, and I want to look at differences in the query sequences (which are much longer than the reference) compared to the reference, and list the results in a table. The position column is related to the position in the reference sequence and not in the query sequence. I want the first AA in the reference sequence to start at 68, not 1.

I find it difficult to add example sequences for this as they tend to be very long, but here is something shorter to work with (not related to the table above):




Expected output:

ID       Reference_AA    Sample_AA   Pos
query1   R               T           70
query2   L               A           72
  • $\begingroup$ What would happen if the reference is found duplicated in the queries? Should they also be returned? Do you want only 1 amino acid substitution? The last part of the first AA of the reference be 68 is easy, just add 67 to any reference position if I understood it correctly. One way would be to use regex with one position changed to ., but then one would need to work how to retrieve the position of the substitution. $\endgroup$ – llrs Jan 25 '18 at 10:28
  • $\begingroup$ If there are multiple differences in a single query, I would want each one to be listed in the output (one for each row). If there is no differences between the query and the reference, nothing needs to be written in the output $\endgroup$ – Haakonkas Jan 25 '18 at 10:30
  • 3
    $\begingroup$ I've provided an accepted answer to this over at stackoverflow $\endgroup$ – heathobrien Jan 25 '18 at 12:44
  • $\begingroup$ I'm voting to close this question as off-topic because it has been cross-posted (and answered) on Stack Overflow. $\endgroup$ – terdon Jan 26 '18 at 13:19

Browse other questions tagged or ask your own question.