1
$\begingroup$

I have a fasta file, like this:

>TRINITY_DN100000_c1_g1::TRINITY_DN100000_c1_g1_i3::g.3039::m.3039 TRINITY_DN100000_c1_g1::TRINITY_DN100000_c1_g1_i3::g.3039  ORF type:complete len:100 (-) TRINITY_DN100000_c1_g1_i3:1027-1326(-)
MVWIKFRGLHRVLTSTPLVKSGKTPSQTWAFLDISVELIVFLFLNVHKSPMPHFKIYSEA
FSEEWSLLWLQYSRHLIQKPKPWQIKIELLHLCCCNRLC*
>TRINITY_DN100000_c1_g6::TRINITY_DN100000_c1_g6_i2::g.84365::m.84365 TRINITY_DN100000_c1_g6::TRINITY_DN100000_c1_g6_i2::g.84365  ORF type:complete len:112 (-) TRINITY_DN100000_c1_g6_i2:379-714(-)
MEMMQEIIPFAREMLSARPSKGTMKVYLVGGTFAVLGIVSGMVEAACSLFPEQEESTLTK
LMEDCLTVTAQNQEPQTFIPEDDEQDAEMEAKAKDLPMFRQRRMSFRAHAS*

if I want to only keep the second header, like this (the amino acid sequence keep unchanged). which command could i USE?

>TRINITY_DN100000_c1_g1_i3
MVWIKFRGLHRVLTSTPLVKSGKTPSQTWAFLDISVELIVFLFLNVHKSPMPHFKIYSEA
FSEEWSLLWLQYSRHLIQKPKPWQIKIELLHLCCCNRLC*
>TRINITY_DN100000_c1_g6_i2
MEMMQEIIPFAREMLSARPSKGTMKVYLVGGTFAVLGIVSGMVEAACSLFPEQEESTLTK
LMEDCLTVTAQNQEPQTFIPEDDEQDAEMEAKAKDLPMFRQRRMSFRAHAS*

I know a command 

> sed 's/::.*//' input > output

can only keep the first header >TRINITY_DN100000_c1_g1 if I want to keep the second header with the isoform information >TRINITY_DN100000_c1_g1_i3, how should I correct this command?

Improve this question
$\endgroup$
1
  • $\begingroup$ You can also do this via eFetch or Biopython directly without regex-ing $\endgroup$
    – M__ ♦
    Oct 7, 2019 at 9:37

2 Answers 2

Reset to default
2
$\begingroup$

Here's one way:

$ sed -E 's/>[^:]*::([^:]*)::.*/>\1/' file.fa 
>TRINITY_DN100000_c1_g1_i3
MVWIKFRGLHRVLTSTPLVKSGKTPSQTWAFLDISVELIVFLFLNVHKSPMPHFKIYSEA
FSEEWSLLWLQYSRHLIQKPKPWQIKIELLHLCCCNRLC*
>TRINITY_DN100000_c1_g6_i2
MEMMQEIIPFAREMLSARPSKGTMKVYLVGGTFAVLGIVSGMVEAACSLFPEQEESTLTK
LMEDCLTVTAQNQEPQTFIPEDDEQDAEMEAKAKDLPMFRQRRMSFRAHAS*

The idea is to match a >, then as many non: as possible ([^:]*), and then capture the string between the first :: and the second ::. Finally, you replace the entire line with the captured string.

Once you're satisfied this does what you need, you can edit the original file with -i:

sed -i -E 's/>[^:]*::([^:]*)::.*/>\1/' file.fa

Alternatively, set awk's field separator to :: and print the second field only on lines that start with >:

$ awk -F'::' '/^>/ ? $0=">"$2 : 1' file.fa 
>TRINITY_DN100000_c1_g1_i3
MVWIKFRGLHRVLTSTPLVKSGKTPSQTWAFLDISVELIVFLFLNVHKSPMPHFKIYSEA
FSEEWSLLWLQYSRHLIQKPKPWQIKIELLHLCCCNRLC*
>TRINITY_DN100000_c1_g6_i2
MEMMQEIIPFAREMLSARPSKGTMKVYLVGGTFAVLGIVSGMVEAACSLFPEQEESTLTK
LMEDCLTVTAQNQEPQTFIPEDDEQDAEMEAKAKDLPMFRQRRMSFRAHAS*
Improve this answer
$\endgroup$
1
$\begingroup$

A 2 step workaround is:

  1. Get rid of everything after the piece of info you want, that is anything starting with :: and goes on with g:

    sed "s/::g.*//" fasta > fasta2

>TRINITY_DN100000_c1_g1::TRINITY_DN100000_c1_g1_i3

  1. Get rid of the "first header" as you call it, after the first step, it is rather straightforward:

    sed "s/>TRINITY_.*::/>/" fasta2 > fasta3

>TRINITY_DN100000_c1_g1_i3

EDIT:

@terdon showed that multiple sed commands can be used when separated by ; in GNU sed (default on Linux):

sed "s/::g.*//; s/>TRINITY_.*::/>/" fasta > modified_fasta
Improve this answer
$\endgroup$
2
  • $\begingroup$ it worked!! thank you! $\endgroup$
    – lili
    Oct 6, 2019 at 22:46
  • 1
    $\begingroup$ There is no need for multiple steps: sed 's/::g.*//; s/>TRINITY_.*::/>/'. Or, if you don't have GNU sed: sed -e 's/::g.*//' -e 's/>TRINITY_.*::/>/'. You could even pipe one sed into another, that would at least avoid the temp files: sed 's/::g.*//' | sed 's/>TRINITY_.*::/>/'. $\endgroup$
    – terdon ♦
    Oct 7, 2019 at 9:30

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.