SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
change order of FASTA seqs, based on ID list ssully General 2 04-23-2013 12:22 AM
Combine two FASTA files Calico Bioinformatics 5 09-26-2012 12:35 PM
how to get specific length sequence from a fasta file entomology Bioinformatics 5 07-12-2012 03:59 PM
fasta sequence: 0 based or 1 based index ardmore Bioinformatics 8 11-15-2011 09:23 AM
Help pls! Picking out specific 'spots' from SRA FASTA files A1_UltiMA Bioinformatics 10 10-13-2010 08:40 AM

Reply
 
Thread Tools
Old 07-08-2013, 07:53 PM   #1
N311V
Member
 
Location: Australia

Join Date: Jul 2013
Posts: 34
Default Combine FASTA files in a specific order based on sequence ID

Hi all,

I frequent these forums often but this is my first post.

I've got a problem that I don't have the scripting skills to solve (nor the time to gain them at the moment).

What I want to do is combine two multi fasta files in a specific order based on the sequence IDs.

For example;

file 1

>seq1
TTTGGATTACAAAGTTATTTAAATCACATGT....
>seq2
GCCGTGCCATTTCAATTACAAATACATAATA....

file 2

>seq1_probe1
CTTTGTCCTTGTCCTTGGTGGCGG....
>seq1_probe2
ATTTCTTCTCATCCTCCTCCTCCTA....
>seq2_probe1
ACTAAAAACTCGTTGAAGAAATCC....
>seq2_probe2
AGGATATAACACACAGCCATCACC....

In need to combined file to look like;

>seq1
TTTGGATTACAAAGTTATTTAAATCACATGT....
>seq1_probe1
CTTTGTCCTTGTCCTTGGTGGCGG....
>seq1
TTTGGATTACAAAGTTATTTAAATCACATGT....
>seq1_probe2
ATTTCTTCTCATCCTCCTCCTCCTA....
>seq2
GCCGTGCCATTTCAATTACAAATACATAATA....
>seq2_probe1
ACTAAAAACTCGTTGAAGAAATCC....
>seq2
GCCGTGCCATTTCAATTACAAATACATAATA....
>seq2_probe2
AGGATATAACACACAGCCATCACC....

Note that only part of file 2's sequence IDs are common to file 1's.

I'd prefer to use perl as that is the language I'm learning but any solution will suffice.

Thanks for reading.
N311V is offline   Reply With Quote
Old 07-08-2013, 08:54 PM   #2
Heisman
Senior Member
 
Location: St. Louis

Join Date: Dec 2010
Posts: 535
Default

I am on my phone and can't type anything elegant (and I don't know perl), but if you want to get the job done with basic linux tools you can look up how to print out every other line in a file with sed (google sed one liners if you can't find it easily), make these separate files, then you can use the paste command followed by the tr command to convert the tabs to new line characters and get what you want. It is ugly but you should be able to figure it out quickly. Use the lines you posted above as test files so you don't waste time practicing with large files.
Heisman is offline   Reply With Quote
Old 07-08-2013, 09:41 PM   #3
Heisman
Senior Member
 
Location: St. Louis

Join Date: Dec 2010
Posts: 535
Default

Here's what I had in mind. Save this in a script, give yourself permission to execute it, and then run it as: ./script file1 file2 output

Code:
#! /bin/bash

file_1=$1
file_2=$2
output=$3

sed -n '1,${p;n}' $file_1 > temp1
sed -n '1,${n;p}' $file_1 > temp2
sed -n '1,${p;n;n;n}' $file_2 > temp3
sed -n '1,${n;p;n;n}' $file_2 > temp4
sed -n '1,${n;n;p;n}' $file_2 > temp5
sed -n '1,${n;n;n;p}' $file_2 > temp6
paste temp1 temp2 temp3 temp4 temp1 temp2 temp5 temp6 | tr '\t' '\n' > $output

rm temp1 temp2 temp3 temp4 temp5 temp6
This is quite inefficient with large files but should introduce some basic commands. You can make it a lot faster by running all of the sed commands together and then having it wait for them to complete prior to putting them together:

Code:
#! /bin/bash

file_1=$1
file_2=$2
output=$3

sed -n '1,${p;n}' $file_1 > temp1 &
pid1=$!
sed -n '1,${n;p}' $file_1 > temp2 &
pid2=$!
sed -n '1,${p;n;n;n}' $file_2 > temp3 &
pid3=$!
sed -n '1,${n;p;n;n}' $file_2 > temp4 &
pid4=$!
sed -n '1,${n;n;p;n}' $file_2 > temp5 &
pid5=$!
sed -n '1,${n;n;n;p}' $file_2 > temp6 &
pid6=$!

wait $pid1 $pid2 $pid3 $pid4 $pid5 $pid6

paste temp1 temp2 temp3 temp4 temp1 temp2 temp5 temp6 | tr '\t' '\n' > $output

rm temp1 temp2 temp3 temp4 temp5 temp6
But obviously with perl you can read in both files and just output the lines in the order you desire. So definitely figure that out too. But it is nice to be able to get stuff done with linux commands while learning how to do things in a much better fashion with a scripting language, so if you can understand how this works that would also be useful.

Last edited by Heisman; 07-08-2013 at 09:44 PM.
Heisman is offline   Reply With Quote
Old 07-09-2013, 12:21 PM   #4
martinghunt
Junior Member
 
Location: Cambridge

Join Date: Jul 2013
Posts: 5
Default

Assuming your files are called 1.fa and 2.fa, this hack will work:

Code:
samtools faidx 2.fa
awk '{id=substr($1,2); getline; for (i=1;i<3;i++){print ">"id; print; system("samtools faidx 2.fa "id"_probe"i)}}' 1.fa
awk is pretty powerful for this kind of thing.

Last edited by martinghunt; 07-09-2013 at 12:23 PM. Reason: didn't need the samtools faidx 1.fa command
martinghunt is offline   Reply With Quote
Old 07-10-2013, 06:15 AM   #5
HMorrison
Senior Member
 
Location: Massachusetts

Join Date: May 2009
Posts: 116
Default

as a one-off solution:

sed -e '$!N;s/\n/\t/' file1 > col1
sed -e '$!N;s/\n/\t/' file2 | sed -e '$!N;s/\n/\t/' > col2
paste col1 col2 | fmt -5

>seq1
TTTGGATTACAAAGTTATTTAAATCACATGT....
>seq1_probe1
CTTTGTCCTTGTCCTTGGTGGCGG....
>seq1_probe2
ATTTCTTCTCATCCTCCTCCTCCTA....
>seq2
GCCGTGCCATTTCAATTACAAATACATAATA....
>seq2_probe1
ACTAAAAACTCGTTGAAGAAATCC....
>seq2_probe2
AGGATATAACACACAGCCATCACC....
HMorrison is offline   Reply With Quote
Old 09-08-2014, 10:45 AM   #6
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

how i can convert
>1...>2....>3...>10000 to >1

and
>1..>2..>3....>10000 for b.fasta to >2 and
same for for all my 5 samples
huma Asif is offline   Reply With Quote
Old 09-08-2014, 10:52 AM   #7
HMorrison
Senior Member
 
Location: Massachusetts

Join Date: May 2009
Posts: 116
Default

Quote:
Originally Posted by huma Asif View Post
how i can convert
>1...>2....>3...>10000 to >1

and
>1..>2..>3....>10000 for b.fasta to >2 and
same for for all my 5 samples
I do not understand the question. Can you explain further?
HMorrison is offline   Reply With Quote
Old 09-08-2014, 11:05 AM   #8
westerman
Rick Westerman
 
Location: Purdue University, Indiana, USA

Join Date: Jun 2008
Posts: 1,104
Default

I agree with @HMorrison -- the question needs to be stated better. That said 'fastx_renamer' will rename FastA files.
westerman is offline   Reply With Quote
Old 09-08-2014, 11:07 AM   #9
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 7,076
Default

This is the parent thread with "some" additional information: http://seqanswers.com/forums/showthread.php?t=46474
GenoMax is offline   Reply With Quote
Old 09-08-2014, 11:19 AM   #10
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

i created fasta from vcf file using target intervals so now in fasta file i have the same number of header as the coordinates in bed

so what i am doing is i want to cat all these sequences
huma Asif is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 04:59 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO