Answers.Body: <p>The <code>exonerate</code> tool comes with some very handy utilities for sequence manipulation. Here, we are interested in <code>fastasubseq</code> and <code>fasterevcomp</code>. First, install <code>exonerate</code> as described <a href="https://www.ebi.ac.uk/about/vertebrate-genomics/software/exonerate" rel="nofollow noreferrer">here</a>. If you're using Debian or Ubuntu or one of their derivatives, it might be enough to run</p>
  
  <pre><code>sudo apt install exonerate
  </code></pre>
  
  <p>Then, follow these steps:</p>
  
  <ol>
  <li><p>Create one fasta file for each sequence in <code>genome.fasta</code></p>
  
  <pre><code>fastaexplode genome.fasta 
  </code></pre></li>
  <li><p>A little preprocessing</p>
  
  <pre><code><span class="math-container">$ awk -F[:-] '{ 
                  if($2&gt;$3){
                      start=$3-1; 
                      len=$2-start" -"
                  }
                  else{
                      start=$2-1; 
                      len=$3-start
                  } 
                  print $</span>1,start,len}' id.txt 
  Contig3 7 8 -
  Contig2 5 5
  Contig1 3 9 -
  </code></pre>
  
  <p>So, that will print the target contig's name, the start coordinate of the target sequence, the length of the target sequence and a <code>-</code> if it is on the reverse strand. The <code>fastasubseq</code> program I will use below starts counting from 0 and not from 1, so we also need to subtract 1 from the coordinates to get them right.</p></li>
  <li><p>Now, feed that into <code>fastafetch</code> and <code>fastarevcomp</code></p>
  
  <pre><code>tmpFile=<span class="math-container">$(mktemp); 
      &gt; subseqs.fa
      awk -F[:-] '{ 
                  if($2&gt;$3){
                      start=$3-1; 
                      len=$2-start" -"
                  }
                  else{
                      start=$2-1; 
                      len=$3-start
                  } 
                  print $1,start,len}' id.txt |
      while read cont start len rev; do 
          fastasubseq "$</span>cont".fa <span class="math-container">$start $</span>len &gt; <span class="math-container">$tmpFile; 
          if [[ -n $</span>rev ]]; then 
              fastarevcomp "<span class="math-container">$tmpFile" &gt;&gt; subseqs.fa; 
          else 
              cat "$</span>tmpFile" &gt;&gt; subseqs.fa; 
      fi
      done
  </code></pre>
  
  <p>This will produce the following file when run on your data:</p>
  
  <pre><code>$ cat subseqs.fa 
  &gt;Contig3:subseq(6,9) [revcomp]
  GATCAGCGC
  &gt;Contig2:subseq(4,6)
  AGCGCT
  &gt;Contig1:subseq(2,10) [revcomp]
  CTTGATGGCC
  </code></pre></li>
  <li><p>Finally, remove the per-contig fasta files</p>
  
  <pre><code>rm Contig*fa
  </code></pre></li>
  </ol>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T06:29:02.870Z
Answers.Text: Thank you so much terdon for your time and help. Your script works awesome. However, I did a mistake in my question. Instead of mentioning Contig3:15-7, Contig1:12-3, I have mentioned Contig3:15-7 and Contig2:5-10.              
                                                                                                                                         I need to reverse complement the sequences, only if it in reverse order/inverse order. Which ever sequences are not in reverse order (Properly order sequences), only to be extracted as its not to be reverse complemented.
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T06:31:16.100Z
Answers.Text: Please accept my sincere apologies for the same. Could you please do the needful.
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T08:21:04.610Z
Answers.Text: @K.Dineshkumar I'm sorry, but I don't understand the problem. The script will only reverse complement if the start position is smaller than the end position. Isn't that what you want? In the example above, `Contig3:15-7` and `Contig1:12-3` have been reverse complimented, but `Contig2:5-10` has not.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T08:39:44.767Z
Answers.Text: @K.Dineshkumar please edit your question and add this. Show us the expected output. It's very hard to read in the comments.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T08:49:05.093Z
Answers.Text: @K.Dineshkumar 
   Never mind, I found the problem. The `fastasubseq` tool uses 0-based coordinates (the first position of a sequence is 0), but you have 1-based (the first position of a sequence is 1). Try the updated answer.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T09:32:43.913Z
Answers.Text: @terdon Please see my edited question, wherein I have included my expected outcome.
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T09:38:14.107Z
Answers.Text: @terdon, The revised/updated script printing one base less in the end of the sequence. For example,  `Contig2:5-10` says, print the sequence from `5 to 10` But the script print only `5-9`.
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T09:56:44.997Z
Answers.Text: @K.Dineshkumar argh! Yes, of course, sorry. Fixed now, try the (last, I hope) update. And thanks for adding the output to your question, now I can check :)
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T11:19:14.567Z
Answers.Text: @terdon, Now it works perfect. I really appreciate your help. Thank you.
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T13:41:39.000Z
Answers.Text: @K.Dineshkumar you're welcome. If this answers solves your issue, please take a moment to [accept it](//$$SITEURL/help/accepted-answer) by clicking on [the checkmark on the left](https://meta.stackexchange.com/a/5235/203101). That will mark the question as answered and is the way that thanks are conveyed on the Stack Exchange sites. You don't need to do it now, you might want to wait in case anyone can give a better solution, but please remember to accept an answer since that marks the question as answered.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-15T17:40:43.953Z
Answers.Text: @K.Dineshkumar thanks for accepting! Feel free to accept one of the other answers though, I won't take it the wrong way :) They give more efficient solutions.
Answers.UserDisplayName: terdon
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-13T15:00:15.163Z
Answers.LastActivityDate: 2019-06-14T09:56:17.523Z
Answers.LastEditDate: 2019-06-14T09:56:17.523Z
Answers.UserDisplayName: terdon
Answers.Body: <p>Since you tagged python, here is a solution that solves this in a single script using Biopython:</p>
  
  <pre><code>import sys
  from Bio import SeqIO
  
  d = {}
  with open('id.txt') as f:
      for line in f:
          contig, coords = line.split(':')
          start, end = map(int, coords.split('-'))
          d[contig] = start, end
  
  for record in SeqIO.parse('genome.fasta', 'fasta'):
      start, end = d[record.id]
      if start &gt; end:
          record.seq = record.seq[end - 1: start].reverse_complement()
      else:
          record.seq = record.seq[start - 1: end]
      SeqIO.write(record, sys.stdout, 'fasta')
  </code></pre>
  
  <p>Output on your example:</p>
  
  <pre><code>&gt;Contig1
  CTTGATGGCC
  &gt;Contig2
  AGCGCT
  &gt;Contig3
  GATCAGCGC
  </code></pre>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-15T18:03:11.653Z
Answers.Text: Thank you @Chris_rands
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-07-06T13:55:24.417Z
Answers.Text: I try to use it for a bacterial genome. But it shows error as follows,          
   `Traceback (most recent call last):
    File "/home/ga/Documents/commands/scripts/amphora/reverse_comple_seq_extrac.py", line 12, in <module>
      start, end = d[record.id]
  KeyError: '1'`
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-07-08T07:37:58.193Z
Answers.Text: @K.Dineshkumar try using `d[record.description]` instead, i can't know what will work without seeing your data though
Answers.UserDisplayName: Chris_Rands
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-07-08T08:40:00.920Z
Answers.Text: My genome of interest can be obtained from the following link `https://www.ncbi.nlm.nih.gov/Traces/wgs/CAHO01?display=contigs`
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-07-08T08:40:26.903Z
Answers.Text: My coordinate is as follows,                                                                                                  `CAHO01000012.1:92353-95295
  CAHO01000017.1:27362-27949
  CAHO01000014.1:11449-12012
  CAHO01000018.1:51156-52928
  CAHO01000030.1:58079-56889
  CAHO01000001.1:60348-62039
  CAHO01000026.1:45206-44508
  `
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-07-08T08:57:45.717Z
Answers.Text: When I tried with `d[record.description]` instead `of d[record.id]`, I end up the error as follows,                                                                                                                                      `Traceback (most recent call last):
    File "/home/ga/Documents/commands/scripts/amphora/reverse_comple_seq_extrac.py", line 12, in <module>
      start, end = d[record.description]
  KeyError: 'CAHO01000001.1 Xanthomonas citri pv. mangiferaeindicae LMG 941 WGS project CAHO00000000 data, contig XMINcontig00001, whole genome shotgun sequence'`
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T14:05:20.187Z
Answers.LastActivityDate: 2019-06-14T14:05:20.187Z
Answers.UserDisplayName: Chris_Rands
Answers.Body: <p>If you're okay with using a Perl script, you could convert your headers into stranded BED6 and use the following script <code>bed2faidxsta.pl</code> with <code>samtools</code>-indexed FASTA to get back stranded sequence:</p>
  
  <p><a href="https://gist.github.com/alexpreynolds/fa9b0f90e181e3b4f640" rel="nofollow noreferrer">https://gist.github.com/alexpreynolds/fa9b0f90e181e3b4f640</a></p>
  
  <p>This script looks at the sixth column of your BED file of chromosome, start, and stop positions, to decide whether to do a reverse complement of the queried sequence. </p>
  
  <p>You could use <code>awk</code> to quickly build such a BED6 file:</p>
  
  <pre><code><span class="math-container">$ awk -vFS='[:-]' '($</span>0~/^&gt;/){ if (<span class="math-container">$2 &lt; $</span>3) { strand = "+"; } else if (<span class="math-container">$2 &gt; $</span>3) { strand = "-"; } else { printf("Coordinates must differ\n"); exit 1; } printf("%s\t%d\t%d\t.\t0\t%c\n", <span class="math-container">$1, $</span>2, $3, strand; }' in.fa &gt; out.bed6
  </code></pre>
  
  <p>To set up indexed FASTA, here's a generic shell-based approach, which uses [Kent utilities][1] <code>fetchChromSizes</code> to get a list of chromosomes for the build <code>mm9</code>:</p>
  
  <pre><code><span class="math-container">$ mkdir mm9
  $</span> cd mm9
  <span class="math-container">$ while read chr; do echo "$</span>chr"; wget -qO- "http://hgdownload.cse.ucsc.edu/goldenpath/mm9/chromosomes/<span class="math-container">$chr.fa.gz" | gunzip -c - &gt; $</span>chr.fa; samtools faidx $chr.fa; done &lt; &lt;(fetchChromSizes mm9 | cut -f1)
  </code></pre>
  
  <p>This can be adjusted for <code>mm10</code>, <code>hg19</code>, <code>hg38</code>, and so on.</p>
  
  <p>Once you have indexed FASTA and the above Perl script, you could do something as easy as the following:</p>
  
  <pre><code>$ ./bed2faidxsta.pl --fastaDir="/path/to/indexed/mm9" &lt; in.bed &gt; out.fa
  </code></pre>
  
  <p>If you're working with a reference genome, the nice thing about this approach is that <code>awk</code>, Perl, and indexed lookups are fast. You can reuse indexed FASTA for other projects that also require sequence lookups by genomic interval.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-15T18:03:29.453Z
Answers.Text: Thank you @Alex Reynolds
Answers.UserDisplayName: Kumar
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2019-06-14T18:41:01.200Z
Answers.LastActivityDate: 2019-06-14T18:41:01.200Z
Answers.UserDisplayName: Alex Reynolds
Body: <p>I have the series of coordinates in <code>id.txt</code> file, whose coordinates sequences are in <code>genome.fasta</code> file.
  The coordinates of id.txt file is shown below,</p>
  
  <pre><code>Contig3:15-7
  Contig2:5-10
  Contig1:12-3
  </code></pre>
  
  <p>The genome.fasta file is shown below,</p>
  
  <pre><code>&gt;Contig1
  AAGGCCATCAAGGACGTGGATGAGGTCGTCAAG
  &gt;Contig2
  ACGAAGCGCTCGCCAAGGCCGAAGAAGAAGGCC 
  &gt;Contig3
  GCTGCGGCGCTGATCCTGGCGGCCCGCGCCGAG
  </code></pre>
  
  <p>I have used the following codes to extract the sequences from <code>genome.fasta</code> based on the coordinates in <code>id.txt</code> by using the following command <code>xargs samtools faidx genome.fasta &lt;id.txt &gt; result.fasta</code>
  It has extracted only <code>Contig2:5-10</code>sequence. Because, the coordinates of <code>Contig2:5-10</code>is in proper order. But <code>Contig3:15-7</code> and <code>Contig1:12-3</code> coordinates are in inverse/reverse order, So samtools could not fetch those sequences.</p>
  
  <p>I need to extract those sequences from <code>Contig1:12-3</code> and <code>Contig3:15-7</code> coordinates and also I want to <code>reverse complement</code> them. I have  larger fasta file and plenty of coordinates sequences to be extracted, So I need to automate this process. So, please help me to automate this process.</p>
  
  <p>The expected output is</p>
  
  <pre><code>&gt;Contig3    
  GATCAGCGC
  &gt;Contig2
  AGCGCT
  &gt;Contig1
  CTTGATGGCC
  </code></pre>
  
  <p>Thank you in advance.</p>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2019-06-13T13:12:25.550Z
Comments.Text: Welcome to the site. Do you have any preference about the tool to be used? What else have you tried to set the coordinates in the right order? Usually this kind of problems are solved by working on every tiny problem one after the other.
Comments.UserDisplayName: llrs
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2019-06-13T13:15:17.153Z
Comments.Text: I have used linux command the following command to keep the coordinates in right order,                                                                                                                                           "awk '$2 > $3 {printf( "%s\t%s\t%s\n", $1, $3, $2); next} id.csv".
Comments.UserDisplayName: Kumar
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2019-06-13T13:17:20.050Z
Comments.Text: So now you have the right order and the positions. What is stopping you from extracting the sequences?
Comments.UserDisplayName: llrs
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2019-06-13T13:18:20.610Z
Comments.Text: By using the above mentioned command I formatted the coordinates in right order and made it compatible for samtools. But still, one more thing I have to do is reverse complement of those right ordered sequences. I have jumbled coordinates, so the command changes the order wherever it needed.
Comments.UserDisplayName: Kumar
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2019-06-13T13:22:47.757Z
Comments.Text: I need to extract the sequence based on the coordinates specified above (even though it is mixed up with inverse ordered and right ordered). But, Now I need to reverse complement the right ordered (right ordered by using above command) sequences. So that I could go for my further analysis. Moreover, I need to automate this process. The reverse complement implementation to be taken place, wherever the inverse coordinates are there.
Comments.UserDisplayName: Kumar
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2019-06-13T14:30:18.180Z
Comments.Text: What do you know biopython? It was designed to make reverse-complement sequences and translate if needed
Comments.UserDisplayName: llrs
Tags: <python><fasta><samtools><shell>
Title: How to reverse complement the DNA sequences for given inverse/reverse coordinates?
ThreadText: # How to reverse complement the DNA sequences for given inverse/reverse coordinates?
  
  (Asked by: username_0 on 2019-06-13)
  
  I have the series of coordinates in `id.txt` file, whose coordinates sequences are in `genome.fasta` file.
  The coordinates of id.txt file is shown below,
  
  ```
  Contig3:15-7
  Contig2:5-10
  Contig1:12-3
  
  ```
  
  The genome.fasta file is shown below,
  
  ```
  >Contig1
  AAGGCCATCAAGGACGTGGATGAGGTCGTCAAG
  >Contig2
  ACGAAGCGCTCGCCAAGGCCGAAGAAGAAGGCC 
  >Contig3
  GCTGCGGCGCTGATCCTGGCGGCCCGCGCCGAG
  
  ```
  
  I have used the following codes to extract the sequences from `genome.fasta` based on the coordinates in `id.txt` by using the following command `xargs samtools faidx genome.fasta <id.txt > result.fasta`
  It has extracted only `Contig2:5-10`sequence. Because, the coordinates of `Contig2:5-10`is in proper order. But `Contig3:15-7` and `Contig1:12-3` coordinates are in inverse/reverse order, So samtools could not fetch those sequences.
  
  I need to extract those sequences from `Contig1:12-3` and `Contig3:15-7` coordinates and also I want to `reverse complement` them. I have larger fasta file and plenty of coordinates sequences to be extracted, So I need to automate this process. So, please help me to automate this process.
  
  The expected output is
  
  ```
  >Contig3    
  GATCAGCGC
  >Contig2
  AGCGCT
  >Contig1
  CTTGATGGCC
  
  ```
  
  Thank you in advance.
  
  **Tags:** python, fasta, samtools, shell
  
  > **Comments:**
  > **username_1**: Welcome to the site. Do you have any preference about the tool to be used? What else have you tried to set the coordinates in the right order? Usually this kind of problems are solved by working on every tiny problem one after the other.  
  > **username_0**: I have used linux command the following command to keep the coordinates in right order,                                                                                                                                           "awk '$2 > $3 {printf( "%s\t%s\t%s\n", $1, $3, $2); next} id.csv".  
  > **username_1**: So now you have the right order and the positions. What is stopping you from extracting the sequences?  
  > **username_0**: By using the above mentioned command I formatted the coordinates in right order and made it compatible for samtools. But still, one more thing I have to do is reverse complement of those right ordered sequences. I have jumbled coordinates, so the command changes the order wherever it needed.  
  > **username_0**: I need to extract the sequence based on the coordinates specified above (even though it is mixed up with inverse ordered and right ordered). But, Now I need to reverse complement the right ordered (right ordered by using above command) sequences. So that I could go for my further analysis. Moreover, I need to automate this process. The reverse complement implementation to be taken place, wherever the inverse coordinates are there.  
  > **username_1**: What do you know biopython? It was designed to make reverse-complement sequences and translate if needed  
  
  ---
  
  ## 3 answers
  
  ### Answer by username_2 on 2019-06-13. Score: 1 | Accepted answer
  
  The `exonerate` tool comes with some very handy utilities for sequence manipulation. Here, we are interested in `fastasubseq` and `fasterevcomp`. First, install `exonerate` as described [here](https://www.ebi.ac.uk/about/vertebrate-genomics/software/exonerate). If you're using Debian or Ubuntu or one of their derivatives, it might be enough to run
  
  ```
  sudo apt install exonerate
  
  ```
  
  Then, follow these steps:
  
  1. Create one fasta file for each sequence in `genome.fasta`
  
     ```
     fastaexplode genome.fasta 
  
     ```
  2. A little preprocessing
  
     ```
     $ awk -F[:-] '{ 
                     if($2>$3){
                         start=$3-1; 
                         len=$2-start" -"
                     }
                     else{
                         start=$2-1; 
                         len=$3-start
                     } 
                     print $1,start,len}' id.txt 
     Contig3 7 8 -
     Contig2 5 5
     Contig1 3 9 -
  
     ```
  
     So, that will print the target contig's name, the start coordinate of the target sequence, the length of the target sequence and a `-` if it is on the reverse strand. The `fastasubseq` program I will use below starts counting from 0 and not from 1, so we also need to subtract 1 from the coordinates to get them right.
  3. Now, feed that into `fastafetch` and `fastarevcomp`
  
     ```
     tmpFile=$(mktemp); 
         > subseqs.fa
         awk -F[:-] '{ 
                     if($2>$3){
                         start=$3-1; 
                         len=$2-start" -"
                     }
                     else{
                         start=$2-1; 
                         len=$3-start
                     } 
                     print $1,start,len}' id.txt |
         while read cont start len rev; do 
             fastasubseq "$cont".fa $start $len > $tmpFile; 
             if [[ -n $rev ]]; then 
                 fastarevcomp "$tmpFile" >> subseqs.fa; 
             else 
                 cat "$tmpFile" >> subseqs.fa; 
         fi
         done
  
     ```
  
     This will produce the following file when run on your data:
  
     ```
     $ cat subseqs.fa 
     >Contig3:subseq(6,9) [revcomp]
     GATCAGCGC
     >Contig2:subseq(4,6)
     AGCGCT
     >Contig1:subseq(2,10) [revcomp]
     CTTGATGGCC
  
     ```
  4. Finally, remove the per-contig fasta files
  
     ```
     rm Contig*fa
  
     ```
  
  > **Comments:**
  > **username_0**: Thank you so much username_2 for your time and help. Your script works awesome. However, I did a mistake in my question. Instead of mentioning Contig3:15-7, Contig1:12-3, I have mentioned Contig3:15-7 and Contig2:5-10.                                                                                                                                                      I need to reverse complement the sequences, only if it in reverse order/inverse order. Which ever sequences are not in reverse order (Properly order sequences), only to be extracted as its not to be reverse complemented.  
  > **username_0**: Please accept my sincere apologies for the same. Could you please do the needful.  
  > **username_2**: @K.Dineshkumar I'm sorry, but I don't understand the problem. The script will only reverse complement if the start position is smaller than the end position. Isn't that what you want? In the example above, `Contig3:15-7` and `Contig1:12-3` have been reverse complimented, but `Contig2:5-10` has not.  
  > **username_2**: @K.Dineshkumar please edit your question and add this. Show us the expected output. It's very hard to read in the comments.  
  > **username_2**: @K.Dineshkumar   Never mind, I found the problem. The `fastasubseq` tool uses 0-based coordinates (the first position of a sequence is 0), but you have 1-based (the first position of a sequence is 1). Try the updated answer.  
  > **username_0**: @username_2 Please see my edited question, wherein I have included my expected outcome.  
  > **username_0**: @username_2, The revised/updated script printing one base less in the end of the sequence. For example,  `Contig2:5-10` says, print the sequence from `5 to 10` But the script print only `5-9`.  
  > **username_2**: @K.Dineshkumar argh! Yes, of course, sorry. Fixed now, try the (last, I hope) update. And thanks for adding the output to your question, now I can check :)  
  > **username_0**: @username_2, Now it works perfect. I really appreciate your help. Thank you.  
  > **username_2**: @K.Dineshkumar you're welcome. If this answers solves your issue, please take a moment to [accept it](//$$SITEURL/help/accepted-answer) by clicking on [the checkmark on the left](https://meta.stackexchange.com/a/5235/203101). That will mark the question as answered and is the way that thanks are conveyed on the Stack Exchange sites. You don't need to do it now, you might want to wait in case anyone can give a better solution, but please remember to accept an answer since that marks the question as answered.  
  > **username_2**: @K.Dineshkumar thanks for accepting! Feel free to accept one of the other answers though, I won't take it the wrong way :) They give more efficient solutions.  
  
  ---
  
  ### Answer by username_3 on 2019-06-14. Score: 2
  
  Since you tagged python, here is a solution that solves this in a single script using Biopython:
  
  ```
  import sys
  from Bio import SeqIO
  
  d = {}
  with open('id.txt') as f:
      for line in f:
          contig, coords = line.split(':')
          start, end = map(int, coords.split('-'))
          d[contig] = start, end
  
  for record in SeqIO.parse('genome.fasta', 'fasta'):
      start, end = d[record.id]
      if start > end:
          record.seq = record.seq[end - 1: start].reverse_complement()
      else:
          record.seq = record.seq[start - 1: end]
      SeqIO.write(record, sys.stdout, 'fasta')
  
  ```
  
  Output on your example:
  
  ```
  >Contig1
  CTTGATGGCC
  >Contig2
  AGCGCT
  >Contig3
  GATCAGCGC
  
  ```
  
  > **Comments:**
  > **username_0**: Thank you @Chris_rands  
  > **username_0**: I try to use it for a bacterial genome. But it shows error as follows,            `Traceback (most recent call last):   File "/home/ga/Documents/commands/scripts/amphora/reverse_comple_seq_extrac.py", line 12, in <module>     start, end = d[record.id] KeyError: '1'`  
  > **username_3**: @K.Dineshkumar try using `d[record.description]` instead, i can't know what will work without seeing your data though  
  > **username_0**: My genome of interest can be obtained from the following link `https://www.ncbi.nlm.nih.gov/Traces/wgs/CAHO01?display=contigs`  
  > **username_0**: My coordinate is as follows,                                                                                                  `CAHO01000012.1:92353-95295 CAHO01000017.1:27362-27949 CAHO01000014.1:11449-12012 CAHO01000018.1:51156-52928 CAHO01000030.1:58079-56889 CAHO01000001.1:60348-62039 CAHO01000026.1:45206-44508 `  
  > **username_0**: When I tried with `d[record.description]` instead `of d[record.id]`, I end up the error as follows,                                                                                                                                      `Traceback (most recent call last):   File "/home/ga/Documents/commands/scripts/amphora/reverse_comple_seq_extrac.py", line 12, in <module>     start, end = d[record.description] KeyError: 'CAHO01000001.1 Xanthomonas citri pv. mangiferaeindicae LMG 941 WGS project CAHO00000000 data, contig XMINcontig00001, whole genome shotgun sequence'`  
  
  ---
  
  ### Answer by username_4 on 2019-06-14. Score: 1
  
  If you're okay with using a Perl script, you could convert your headers into stranded BED6 and use the following script `bed2faidxsta.pl` with `samtools`-indexed FASTA to get back stranded sequence:
  
  <https://gist.github.com/alexpreynolds/fa9b0f90e181e3b4f640>
  
  This script looks at the sixth column of your BED file of chromosome, start, and stop positions, to decide whether to do a reverse complement of the queried sequence.
  
  You could use `awk` to quickly build such a BED6 file:
  
  ```
  $ awk -vFS='[:-]' '($0~/^>/){ if ($2 < $3) { strand = "+"; } else if ($2 > $3) { strand = "-"; } else { printf("Coordinates must differ\n"); exit 1; } printf("%s\t%d\t%d\t.\t0\t%c\n", $1, $2, $3, strand; }' in.fa > out.bed6
  
  ```
  
  To set up indexed FASTA, here's a generic shell-based approach, which uses [Kent utilities][1] `fetchChromSizes` to get a list of chromosomes for the build `mm9`:
  
  ```
  $ mkdir mm9
  $ cd mm9
  $ while read chr; do echo "$chr"; wget -qO- "http://hgdownload.cse.ucsc.edu/goldenpath/mm9/chromosomes/$chr.fa.gz" | gunzip -c - > $chr.fa; samtools faidx $chr.fa; done < <(fetchChromSizes mm9 | cut -f1)
  
  ```
  
  This can be adjusted for `mm10`, `hg19`, `hg38`, and so on.
  
  Once you have indexed FASTA and the above Perl script, you could do something as easy as the following:
  
  ```
  $ ./bed2faidxsta.pl --fastaDir="/path/to/indexed/mm9" < in.bed > out.fa
  
  ```
  
  If you're working with a reference genome, the nice thing about this approach is that `awk`, Perl, and indexed lookups are fast. You can reuse indexed FASTA for other projects that also require sequence lookups by genomic interval.
  
  > **Comments:**
  > **username_0**: Thank you @username_4  
  
  ---
  
  