Question: How to append two fasta files?
0
Entering edit mode

I have two fasta files as shown below,

File:1

Contig_1:90600-91187 AAGGCCATCAAGGACGTGGATGAGGTCGTCAAGGGCAAGGAACAGGAATTGATGACGGTC Contig_98:35323-35886 GACGAAGCGCTCGCCAAGGCCGAAGAAGAAGGCCTGGATCTGGTCGAAATCCAGCCGCAG Contig_24:26615-28387 GCTGCGGCGCTGATCCTGGCGGCCCGCGCCGAGGAGATCGCCCGTTTGGAGCGCGGCGAA

File:2

Contig_1:90600-91187
GACCGTCATCAATTCCTGTTCCTTGCCCTTGACGACCTCATCCACGTCCTTGATGGCCTT
Contig_24:26615-28387
TTCGCCGCGCTCCAAACGGGCGATCTCCTCGGCGCGGGCCGCCAGGATCAGCGCCG

Both files are having the same fasta headers but vary in terms of sequences. Hence, I need to replace File:2 sequences in File:1 as shown below,

Expected outcome:

Contig_1:90600-91187
GACCGTCATCAATTCCTGTTCCTTGCCCTTGACGACCTCATCCACGTCCTTGATGGCCTT
Contig_98:35323-35886
GACGAAGCGCTCGCCAAGGCCGAAGAAGAAGGCCTGGATCTGGTCGAAATCCAGCCG
Contig_24:26615-28387
TTCGCCGCGCTCCAAACGGGCGATCTCCTCGGCGCGGGCCGCCAGGATCAGCGCCG

I tried with cat command, but it is concatenating all the sequences instead of replacing the sequences as I mentioned above.

ADD COMMENTlink 7 months ago Dineshkumar K • 10 • updated 7 months ago Corentin • 290
0
Entering edit mode

Hi,

You can use "samtools faidx"

samtools faidx file1 Contig_1:90600-91187  > result.fa

samtools faidx file2 Contig_98:35323-35886 >> result.fa
samtools faidx file2 Contig_24:26615-28387 >> result.fa

However, the name format "contig:start-stop" could cause some issues, faidx might interpret it as a command to extract just a portion of the sequences. You can rename the contigs as "name_start_stop" to avoid this.

ADD COMMENTlink 7 months ago Corentin • 290
Entering edit mode
0

Thank you Corentin. However, I have large fasta file, in that case it would be tedious to implement samtools faidx.

ADD REPLYlink 7 months ago
Dineshkumar K
• 10
Entering edit mode
0

There is an option to give a list of regions to faidx:

-r, --region-file FILE   File of regions.  Format is chr:from-to. One per line.

You could create two files, one with a list of regions from file1, the other from file2 and run two faidx commands

ADD REPLYlink 7 months ago
Corentin
• 290
Entering edit mode
0

Thank you Corentin, Is it possible to extract the sequences based on the reverse order coordinates by samtools faidx. For example, Instead of proper order

Contig_98:35323-35886

If it like this,

Contig_98:35886-35323

Will it work?

ADD REPLYlink 7 months ago
Dineshkumar K
• 10
Entering edit mode
1

I am not exactly sure what you want to do, if you want the reverse complement, faidx also have an option for that:

   -i, --reverse-complement
    Output the sequence as the reverse complement. `When this option is used, “/rc” will be appended to the sequence names. To turn this off or change the string appended, use the --mark-strand option`.

You can see the documentation here : http://www.htslib.org/doc/samtools.html

ADD REPLYlink 7 months ago
Corentin
• 290
Entering edit mode
0

Dear Corentin, If the coordinates are in reverse order, will samtools work? Let me put in this way, I have the series of coordinates to be extracted from the given data,

Contig_98:35323-35886
Contig_1:90600-91187
Contig_24:26615-28387
Contig_91:690-450

The above mentioned coordinates are well formatted except the last one Contig_91:690-450 . Because, the last coordinate is in reverse order. In that case will samtool works?

ADD REPLYlink 7 months ago
Dineshkumar K
• 10
Entering edit mode
1

Apparently faidx is giving this error if the coordinates are reversed:

[faidx] Zero length sequence:
ADD REPLYlink 7 months ago
Corentin
• 290
Entering edit mode
0

Thank you Corentin for your time and help. I have implemented your suggestions and all works well.

ADD REPLYlink 7 months ago
Dineshkumar K
• 10

Login before adding your answer.

Powered by the version 1.8