text: user292: I have issues working with pybedtools due to strand information being lost after extracting and merging transcript coordinates from a gtf file.
  It seems that the solution to preserve strand information depends on the version of bedtools, so I would like to have my code detect bedtools version and adapt its behaviour accordingly.
  On a machine with bedtools v2.27.0:
  ids = ['WBGene00022277', 'WBGene00007068']
  
  def in_ids(itvl):
      return itvl[2] == "transcript" and itvl["gene_id"] in ids
  
  all_genes = BedTool("genes.gtf")
  my_genes = all_genes.filter(in_ids).saveas()
  print(my_genes)
  
  I       ensembl transcript      4119    10230   .       -       .       gene_biotype "protein_coding";gene_id "WBGene00022277";gene_name "homt-1";gene_source "ensembl";gene_version "1";p_id "P18791";transcript_biotype "protein_coding";transcript_id "Y74C9A.3";transcript_name "Y74C9A.3";transcript_source "ensembl";transcript_version "1";tss_id "TSS21161";
  X       ensembl transcript      17714968        17718714        .       +       .       gene_biotype "protein_coding";gene_id "WBGene00007068";gene_name "cTel55X.1";gene_source "ensembl";gene_version "1";p_id "P8395";transcript_biotype "protein_coding";transcript_id "cTel55X.1b";transcript_name "cTel55X.1b";transcript_source "ensembl";transcript_version "1";tss_id "TSS18767";
  X       ensembl transcript      17714968        17718717        .       +       .       gene_biotype "protein_coding";gene_id "WBGene00007068";gene_name "cTel55X.1";gene_source "ensembl";gene_version "1";p_id "P8874";transcript_biotype "protein_coding";transcript_id "cTel55X.1a";transcript_name "cTel55X.1a";transcript_source "ensembl";transcript_version "1";tss_id "TSS18767";
  
  merged = my_genes.sort().merge(s=True).saveas()
  print(merged)
  
  I       4118    10230
  X       17714967        17718717
  
  (Missing strand information)
  merged = my_genes.sort().merge(s=True, c=7, o="distinct").saveas()
  print(merged)
  
  I       4118    10230   -
  X       17714967        17718717        +
  
  On a machine with bedtools v2.26.0:
  ids = ['WBGene00022277', 'WBGene00007068']
  
  def in_ids(itvl):
      return itvl[2] == "transcript" and itvl["gene_id"] in ids
  
  all_genes = BedTool("genes.gtf")
  my_genes = all_genes.filter(in_ids).saveas()
  print(my_genes)
  
  I       ensembl transcript      4119    10230   .       -       .       gene_biotype "protein_coding";gene_id "WBGene00022277";gene_name "homt-1";gene_source "ensembl";gene_version "1";p_id "P18791";transcript_biotype "protein_coding";transcript_id "Y74C9A.3";transcript_name "Y74C9A.3";transcript_
  source "ensembl";transcript_version "1";tss_id "TSS21161";
  X       ensembl transcript      17714968        17718714        .       +       .       gene_biotype "protein_coding";gene_id "WBGene00007068";gene_name "cTel55X.1";gene_source "ensembl";gene_version "1";p_id "P8395";transcript_biotype "protein_coding";transcript_id "cTel55X.1b";transcript_name "c
  Tel55X.1b";transcript_source "ensembl";transcript_version "1";tss_id "TSS18767";
  X       ensembl transcript      17714968        17718717        .       +       .       gene_biotype "protein_coding";gene_id "WBGene00007068";gene_name "cTel55X.1";gene_source "ensembl";gene_version "1";p_id "P8874";transcript_biotype "protein_coding";transcript_id "cTel55X.1a";transcript_name "c
  Tel55X.1a";transcript_source "ensembl";transcript_version "1";tss_id "TSS18767";
  
  merged = my_genes.sort().merge(s=True).saveas()
  print(merged)
  
  I       4118    10230   -
  X       17714967        17718717        +
  
  merged = my_genes.sort().merge(s=True, c=7, o="distinct").saveas()
  print(merged)
  
  I       4118    10230   -       -
  X       17714967        17718717        +       +
  
  (Duplicate strand information)
  Given that I work on linux systems, how can I detect bedtools version from within python?
  Is there a better way to do what I want?
  (Ideally, I would like to have my resulting bed consist in chrom, start, stop, name, score, strand fields, where the score would just be "." and the name would be the content of itvl["gene_id"])
  
  Edit (05/05/2021): Regarding my last comment above, I found out that it is easy to convert gff into bed using pybedtools.featurefuncs.gff2bed (tested with bedtools 2.30):
  from pybedtools import BedTool
  from pybedtools.featurefuncs import gff2bed
  
  ids = ['WBGene00022277', 'WBGene00007068']
  
  def in_ids(itvl):
      return itvl[2] == "transcript" and itvl["gene_id"] in ids
  
  all_genes = BedTool("genes.gtf")
  my_genes = all_genes.filter(in_ids).each(gff2bed, "gene_id").saveas()
  print(my_genes)
  
  I       4118    10230   WBGene00022277  .       -
  X       17714967        17718714        WBGene00007068  .       +
  X       17714967        17718717        WBGene00007068  .       +
  
  merged = my_genes.sort().merge(s=True, c=[4, 5, 6], o="distinct").saveas()
  print(merged)
  
  I       4118    10230   WBGene00022277  .       -
  X       17714967        17718717        WBGene00007068  .       +
  
  This ensures proper bed format.
  However, with bedtools 2.26, s=True makes strand information appear as a 4th column, so it is more complicated to obtain correct column order after merging.
  
  user77: pybedtools assumes that bedtools is in your path and bedtools itself will return the version with bedtools --version. So:
  import subprocess
  subprocess.check_output(['bedtools', '--version'], text=True).strip().split()[1]
  
  For me that returns 'v2.30.0'.
  
original_text: 292: <p>I have issues working with pybedtools due to strand information being lost after extracting and merging transcript coordinates from a gtf file.</p>
  <p>It seems that the solution to preserve strand information depends on the version of bedtools, so I would like to have my code detect bedtools version and adapt its behaviour accordingly.</p>
  <p>On a machine with bedtools v2.27.0:</p>
  <pre class="lang-python prettyprint-override"><code>ids = ['WBGene00022277', 'WBGene00007068']
  
  def in_ids(itvl):
      return itvl[2] == &quot;transcript&quot; and itvl[&quot;gene_id&quot;] in ids
  
  all_genes = BedTool(&quot;genes.gtf&quot;)
  my_genes = all_genes.filter(in_ids).saveas()
  print(my_genes)
  </code></pre>
  <pre><code>I       ensembl transcript      4119    10230   .       -       .       gene_biotype &quot;protein_coding&quot;;gene_id &quot;WBGene00022277&quot;;gene_name &quot;homt-1&quot;;gene_source &quot;ensembl&quot;;gene_version &quot;1&quot;;p_id &quot;P18791&quot;;transcript_biotype &quot;protein_coding&quot;;transcript_id &quot;Y74C9A.3&quot;;transcript_name &quot;Y74C9A.3&quot;;transcript_source &quot;ensembl&quot;;transcript_version &quot;1&quot;;tss_id &quot;TSS21161&quot;;
  X       ensembl transcript      17714968        17718714        .       +       .       gene_biotype &quot;protein_coding&quot;;gene_id &quot;WBGene00007068&quot;;gene_name &quot;cTel55X.1&quot;;gene_source &quot;ensembl&quot;;gene_version &quot;1&quot;;p_id &quot;P8395&quot;;transcript_biotype &quot;protein_coding&quot;;transcript_id &quot;cTel55X.1b&quot;;transcript_name &quot;cTel55X.1b&quot;;transcript_source &quot;ensembl&quot;;transcript_version &quot;1&quot;;tss_id &quot;TSS18767&quot;;
  X       ensembl transcript      17714968        17718717        .       +       .       gene_biotype &quot;protein_coding&quot;;gene_id &quot;WBGene00007068&quot;;gene_name &quot;cTel55X.1&quot;;gene_source &quot;ensembl&quot;;gene_version &quot;1&quot;;p_id &quot;P8874&quot;;transcript_biotype &quot;protein_coding&quot;;transcript_id &quot;cTel55X.1a&quot;;transcript_name &quot;cTel55X.1a&quot;;transcript_source &quot;ensembl&quot;;transcript_version &quot;1&quot;;tss_id &quot;TSS18767&quot;;
  </code></pre>
  <pre class="lang-python prettyprint-override"><code>merged = my_genes.sort().merge(s=True).saveas()
  print(merged)
  </code></pre>
  <pre><code>I       4118    10230
  X       17714967        17718717
  </code></pre>
  <p>(<strong>Missing strand information</strong>)</p>
  <pre class="lang-python prettyprint-override"><code>merged = my_genes.sort().merge(s=True, c=7, o=&quot;distinct&quot;).saveas()
  print(merged)
  </code></pre>
  <pre><code>I       4118    10230   -
  X       17714967        17718717        +
  </code></pre>
  <p>On a machine with bedtools v2.26.0:</p>
  <pre class="lang-python prettyprint-override"><code>ids = ['WBGene00022277', 'WBGene00007068']
  
  def in_ids(itvl):
      return itvl[2] == &quot;transcript&quot; and itvl[&quot;gene_id&quot;] in ids
  
  all_genes = BedTool(&quot;genes.gtf&quot;)
  my_genes = all_genes.filter(in_ids).saveas()
  print(my_genes)
  </code></pre>
  <pre><code>I       ensembl transcript      4119    10230   .       -       .       gene_biotype &quot;protein_coding&quot;;gene_id &quot;WBGene00022277&quot;;gene_name &quot;homt-1&quot;;gene_source &quot;ensembl&quot;;gene_version &quot;1&quot;;p_id &quot;P18791&quot;;transcript_biotype &quot;protein_coding&quot;;transcript_id &quot;Y74C9A.3&quot;;transcript_name &quot;Y74C9A.3&quot;;transcript_
  source &quot;ensembl&quot;;transcript_version &quot;1&quot;;tss_id &quot;TSS21161&quot;;
  X       ensembl transcript      17714968        17718714        .       +       .       gene_biotype &quot;protein_coding&quot;;gene_id &quot;WBGene00007068&quot;;gene_name &quot;cTel55X.1&quot;;gene_source &quot;ensembl&quot;;gene_version &quot;1&quot;;p_id &quot;P8395&quot;;transcript_biotype &quot;protein_coding&quot;;transcript_id &quot;cTel55X.1b&quot;;transcript_name &quot;c
  Tel55X.1b&quot;;transcript_source &quot;ensembl&quot;;transcript_version &quot;1&quot;;tss_id &quot;TSS18767&quot;;
  X       ensembl transcript      17714968        17718717        .       +       .       gene_biotype &quot;protein_coding&quot;;gene_id &quot;WBGene00007068&quot;;gene_name &quot;cTel55X.1&quot;;gene_source &quot;ensembl&quot;;gene_version &quot;1&quot;;p_id &quot;P8874&quot;;transcript_biotype &quot;protein_coding&quot;;transcript_id &quot;cTel55X.1a&quot;;transcript_name &quot;c
  Tel55X.1a&quot;;transcript_source &quot;ensembl&quot;;transcript_version &quot;1&quot;;tss_id &quot;TSS18767&quot;;
  </code></pre>
  <pre class="lang-python prettyprint-override"><code>merged = my_genes.sort().merge(s=True).saveas()
  print(merged)
  </code></pre>
  <pre><code>I       4118    10230   -
  X       17714967        17718717        +
  </code></pre>
  <pre class="lang-python prettyprint-override"><code>merged = my_genes.sort().merge(s=True, c=7, o=&quot;distinct&quot;).saveas()
  print(merged)
  </code></pre>
  <pre><code>I       4118    10230   -       -
  X       17714967        17718717        +       +
  </code></pre>
  <p>(<strong>Duplicate strand information</strong>)</p>
  <p>Given that I work on linux systems, <strong>how can I detect bedtools version from within python</strong>?</p>
  <p><strong>Is there a better way</strong> to do what I want?</p>
  <p>(Ideally, I would like to have my resulting bed consist in chrom, start, stop, name, score, strand fields, where the score would just be &quot;.&quot; and the name would be the content of <code>itvl[&quot;gene_id&quot;]</code>)</p>
  <hr />
  <p>Edit (05/05/2021): Regarding my last comment above, I found out that it is easy to convert gff into bed using <a href="https://daler.github.io/pybedtools/pybedtools.featurefuncs.gff2bed.html#pybedtools.featurefuncs.gff2bed" rel="nofollow noreferrer">pybedtools.featurefuncs.gff2bed</a> (tested with bedtools 2.30):</p>
  <pre class="lang-python prettyprint-override"><code>from pybedtools import BedTool
  from pybedtools.featurefuncs import gff2bed
  
  ids = ['WBGene00022277', 'WBGene00007068']
  
  def in_ids(itvl):
      return itvl[2] == &quot;transcript&quot; and itvl[&quot;gene_id&quot;] in ids
  
  all_genes = BedTool(&quot;genes.gtf&quot;)
  my_genes = all_genes.filter(in_ids).each(gff2bed, &quot;gene_id&quot;).saveas()
  print(my_genes)
  </code></pre>
  <pre><code>I       4118    10230   WBGene00022277  .       -
  X       17714967        17718714        WBGene00007068  .       +
  X       17714967        17718717        WBGene00007068  .       +
  </code></pre>
  <pre class="lang-python prettyprint-override"><code>merged = my_genes.sort().merge(s=True, c=[4, 5, 6], o=&quot;distinct&quot;).saveas()
  print(merged)
  </code></pre>
  <pre><code>I       4118    10230   WBGene00022277  .       -
  X       17714967        17718717        WBGene00007068  .       +
  </code></pre>
  <p>This ensures proper bed format.</p>
  <p>However, with bedtools 2.26, <code>s=True</code> makes strand information appear as a 4th column, so it is more complicated to obtain correct column order after merging.</p>
   77: <p>pybedtools assumes that bedtools is in your path and bedtools itself will return the version with <code>bedtools --version</code>. So:</p>
  <pre><code>import subprocess
  subprocess.check_output(['bedtools', '--version'], text=True).strip().split()[1]
  </code></pre>
  <p>For me that returns <code>'v2.30.0'</code>.</p>
  