<p>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.</p>

<p>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.</p>

<p>I want to end up with something like this:</p>

<pre><code>Isolate ID    Reference_AA   Sample_AA   Pos
1             S              T           254
2             T              D           200
3             L              A           230
</code></pre>

<p>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.</p>

<p>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):</p>

<pre><code>&gt;reference
VGRALPDVR

&gt;query1
KSSYLDYAMSVIVGTALPDVRDGLKPVHRRVLY

&gt;query2
ELKSSYLDYAMSVIVGRAAPDVRDGLKPV

Expected output:

ID       Reference_AA    Sample_AA   Pos
query1   R               T           70
query2   L               A           72
</code></pre>

CC BY-SA 3.0
2018-01-25T10:28:27.650Z
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.
llrs
CC BY-SA 3.0
2018-01-25T10:30:41.490Z
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
Haakonkas
CC BY-SA 3.0
2018-01-25T12:44:53.690Z
I've provided an accepted answer to this over at [stackoverflow](https://stackoverflow.com/questions/48402460/identifying-amino-acid-substitutions-from-local-alignments-in-r)
heathobrien
CC BY-SA 3.0
2018-01-26T13:19:37.450Z
I'm voting to close this question as off-topic because it has been cross-posted (and answered) on [Stack Overflow](https://stackoverflow.com/q/48402460/1081936).
terdon
<r><bioconductor>
Identifying amino acid substitutions from local alignments in R
# Identifying amino acid substitutions from local alignments in R

(Asked by: username_0 on 2018-01-23)

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):
>reference
VGRALPDVR

>query1
KSSYLDYAMSVIVGTALPDVRDGLKPVHRRVLY

>query2
ELKSSYLDYAMSVIVGRAAPDVRDGLKPV

Expected output:

ID       Reference_AA    Sample_AA   Pos
query1   R               T           70
query2   L               A           72

**Tags:** r, bioconductor
> **Comments:**
> **username_1**: 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.  
> **username_0**: 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  
> **username_2**: I've provided an accepted answer to this over at [stackoverflow](https://stackoverflow.com/questions/48402460/identifying-amino-acid-substitutions-from-local-alignments-in-r)  
> **username_3**: I'm voting to close this question as off-topic because it has been cross-posted (and answered) on [Stack Overflow](https://stackoverflow.com/q/48402460/1081936).  

---

## 0 answers

