How to create perl script for comparing two sequence files and print the matching sequences in output file
1
2
Entering edit mode
7.0 years ago
vinij80 ▴ 30

Hi,

I have two sequence file. I want only the matching sequences which is present on both files and save that output in third file. Can you please help me for creating perl program for this?

FILE1:

>Contig1
TTCAAAAACTCATATGGGTGGTACAATGCGTCTTGGATCTAGGAGAACATATTTTCAAGTTGCAGATTGTAAATCTGCAAAATTATATGGTAACCAGAGCTTTGTAGATGAGAGGCATCGACACAGATATGAGGTGAACCCCGACATGGTGCAGC
>Contig2
GACTTGAAGATGCTGGTCTTTCTTTCACTGGCAAAGATGAAAGTGGTCATCGCATGGAGATTGTTGAGCTGCCGAGTCATCCTTACTTCATCGGAGTTCAATTTCATCCAGAATTTAAATCAAGGCCAGGAACCCCTTCAGCCCTGTTT
>Contig3
CTAGGACTTATAGCCGCAGCAACTGGGCAACTTGAAACTCTCTTGAAGAAGGGTGTTCCCAAAACATGGGGGTTGAGCAATGGTACGTCAGGACTAAAATCACATCGATATGTAAATGGGACAAAACTGTTTAATGGATCATTAGATG
>Contig4
GCATTTATTGCAATGGGAATGGTATACATGTTTAAAGGAAACAGTAACATATGTTGTGGGCGCTTGGCCCCGGATTTTTGATAATCAAATTTTGCTACTGCATTTTTTTTAAAG
>Contig5
CCCCCCCTTATTTGTCGTTTTTGATAATCAAATTTTGCTACTGCATTTTTTTTAAAG

FILE2:

>Contig1
TTCAAAAACTCATATGGGTGGTACAATGCGTCTTGGATCTAGGAGAACATATTTTCAAGTTGCAGATTGTAAATCTGCAAAATTATATGGTAACCAGAGCTTTGTAGATGAGAGGCATCGACACAGATATGAGGTGAACCCCGACATGGTGCAGC
>Contig3
CTAGGACTTATAGCCGCAGCAACTGGGCAACTTGAAACTCTCTTGAAGAAGGGTGTTCCCAAAACATGGGGGTTGAGCAATGGTACGTCAGGACTAAAATCACATCGATATGTAAATGGGACAAAACTGTTTAATGGATCATTAGATG
>Contig4
GCATTTATTGCAATGGGAATGGTATACATGTTTAAAGGAAACAGTAACATATGTTGTGGGCGCTTGGCCCCGGATTTTTGATAATCAAATTTTGCTACTGCATTTTTTTTAAAG
>Contig6
GCATTTATTGCAATGTTTTGATAATCAAATTTTGCTACTGCATTTTTTTTAAAGCCAGAGCTTTGTAGATGAGAGGCATGGTACAATGCGTCTTG

EXPECTED OUTPUT:

>Contig1
TTCAAAAACTCATATGGGTGGTACAATGCGTCTTGGATCTAGGAGAACATATTTTCAAGTTGCAGATTGTAAATCTGCAAAATTATATGGTAACCAGAGCTTTGTAGATGAGAGGCATCGACACAGATATGAGGTGAACCCCGACATGGTGCAGC
>Contig3
CTAGGACTTATAGCCGCAGCAACTGGGCAACTTGAAACTCTCTTGAAGAAGGGTGTTCCCAAAACATGGGGGTTGAGCAATGGTACGTCAGGACTAAAATCACATCGATATGTAAATGGGACAAAACTGTTTAATGGATCATTAGATG
>Contig4
GCATTTATTGCAATGGGAATGGTATACATGTTTAAAGGAAACAGTAACATATGTTGTGGGCGCTTGGCCCCGGATTTTTGATAATCAAATTTTGCTACTGCATTTTTTTTAAAG
sequence • 3.3k views
ADD COMMENT
2
Entering edit mode

What do you want to do is something basic. So, if you are learning Perl and bioinformatics you may have a look to this perl guide. You have all information in this guide in order to reach your goal.

Then, if we talk about algorithms, what you have to use is something like this:

  1. Open the two files to read sequences.
  2. For each sequence in FILE 1, you have to compare the sequence header in all sequences in FILE 2.
  3. If headers are the same (or sequence, you can compare both, but I think that if the sequence and headers are the same, it is better to only compare header), you have to write the header and the associated sequence in a file.
  4. Close files, and you have what you need! ;)
ADD REPLY
0
Entering edit mode

I don't know how to write this program. Can you please help me?

ADD REPLY
0
Entering edit mode

Use BioPython or BioRuby or a software like DNAStar

ADD REPLY
0
Entering edit mode

Hello vinij80!

It appears that your post has been cross-posted to another site: cross posted: http://seqanswers.com/forums/showthread.php?t=63159

This is typically not recommended as it runs the risk of annoying people in both communities.

ADD REPLY
2
Entering edit mode
7.0 years ago
bioguy24 ▴ 220

If awk can be used:

awk 'NR==FNR{a[$0];next}$0 in a{print $0}' file1.txt file2.txt
>Contig1
TTCAAAAACTCATATGGGTGGTACAATGCGTCTTGGATCTAGGAGAACATATTTTCAAGTTGCAGATTGTAAATCTGCAAAATTATATGGTAACCAGAGCTTTGTAGATGAGAGGCATCGACACAGATATGAGGTGAACCCCGACATGGTGCAGC
>Contig3
CTAGGACTTATAGCCGCAGCAACTGGGCAACTTGAAACTCTCTTGAAGAAGGGTGTTCCCAAAACATGGGGGTTGAGCAATGGTACGTCAGGACTAAAATCACATCGATATGTAAATGGGACAAAACTGTTTAATGGATCATTAGATG
>Contig4
GCATTTATTGCAATGGGAATGGTATACATGTTTAAAGGAAACAGTAACATATGTTGTGGGCGCTTGGCCCCGGATTTTTGATAATCAAATTTTGCTACTGCATTTTTTTTAAAG

A perl equivalent

perl -ne 'print if ($seen{$_} .= @ARGV) =~ /10$/'  file1.txt file2.txt
>Contig1
TTCAAAAACTCATATGGGTGGTACAATGCGTCTTGGATCTAGGAGAACATATTTTCAAGTTGCAGATTGTAAATCTGCAAAATTATATGGTAACCAGAGCTTTGTAGATGAGAGGCATCGACACAGATATGAGGTGAACCCCGACATGGTGCAGC
>Contig3
CTAGGACTTATAGCCGCAGCAACTGGGCAACTTGAAACTCTCTTGAAGAAGGGTGTTCCCAAAACATGGGGGTTGAGCAATGGTACGTCAGGACTAAAATCACATCGATATGTAAATGGGACAAAACTGTTTAATGGATCATTAGATG
>Contig4
GCATTTATTGCAATGGGAATGGTATACATGTTTAAAGGAAACAGTAACATATGTTGTGGGCGCTTGGCCCCGGATTTTTGATAATCAAATTTTGCTACTGCATTTTTTTTAAAG
ADD COMMENT
0
Entering edit mode

Thankyou ....But i want the output in a separate text file

ADD REPLY
0
Entering edit mode

Just redirect the output to a text file then you will get a seperate file with what I showed below the script

awk 'NR==FNR{a[$0];next}$0 in a{print $0}' file1.txt file2.txt > output.txt
perl -ne 'print if ($seen{$_} .= @ARGV) =~ /10$/'  file1.txt file2.txt > output.txt
ADD REPLY

Login before adding your answer.

Traffic: 1926 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6