Merge pull request #168 from JureZmrzlikar/bedGraph

bedgraph: Add BED to bedGraph converter

iCount/cli.py

@@ -366,6 +366,9 @@ def main():
     make_parser_from_function(
         iCount.analysis.summary.make_summary_report, subparsers)
 
+    # File converters:
+    make_parser_from_function(iCount.files.bedgraph.bed2bedgraph, subparsers)
+
     # Examples:
     make_parser_from_function(
         iCount.examples.run, subparsers)

iCount/files/__init__.py

@@ -13,6 +13,9 @@
 .. automodule:: iCount.files.bed
    :members:
 
+.. automodule:: iCount.files.bedgraph
+   :members:
+
 .. automodule:: iCount.files.fastq
    :members:
 
@@ -48,6 +51,7 @@
 import iCount
 
 from . import bed
+from . import bedgraph
 from . import fasta
 from . import fastq
 

iCount/files/bedgraph.py

@@ -0,0 +1,59 @@
+""".. Line to protect from pydocstyle D205, D400.
+
+Bedgraph conversion
+-------------------
+
+Convert from BED6 to bedGraph format.
+
+"""
+import logging
+
+import pybedtools
+
+import iCount
+
+LOGGER = logging.getLogger(__name__)
+
+
+def bed2bedgraph(bed, bedgraph, name='User Track', description='User Supplied Track'):
+    """
+    Convert from BED6 to bedGraph format.
+
+    For further explanation of ``name`` and ``description`` parameters
+    (and also others, that are not yet supported) see:
+    https://genome.ucsc.edu/goldenPath/help/customTrack.html#TRACK
+
+    Parameters
+    ----------
+    bed : str
+        Input BED6 file.
+    bedgraph : str
+        Output bedGraph file.
+    name : str
+        Track label. Can consist of up to 15 characters.
+    description : str
+        Track description. Can consist of up to 60 characters.
+
+    Returns
+    -------
+    None
+
+    """
+    iCount.log_inputs(LOGGER, level=logging.INFO)  # pylint: disable=protected-access
+    LOGGER.info('Checking input parameters...')
+    assert bed.endswith(('.bed', '.bed.gz'))
+
+    header = [
+        'track',
+        'type=bedGraph',
+        'name="{}"'.format(name[:15]),
+        'description="{}"'.format(description[:60]),
+    ]
+
+    LOGGER.info('Converting %s to %s', bed, bedgraph)
+    with open(bedgraph, 'wt') as outfile:
+        outfile.write(' '.join(map(str, header)) + '\n')
+        for line in pybedtools.BedTool(bed):
+            bg_line = [line.chrom, line.start, line.stop, '{}{}'.format(line.strand, line.score)]
+            outfile.write('\t'.join(map(str, bg_line)) + '\n')
+    LOGGER.info('Done.')

iCount/tests/test_cli.py

@@ -26,7 +26,7 @@ def setUp(self):
             ['1', '16', '17', '.', '5', '+'],
             ['1', '14', '15', '.', '5', '+'],
             ['1', '15', '16', '.', '5', '+'],
-        ])
+        ], extension='bed')
 
         self.peaks = make_file_from_list([
             ['1', '15', '16', '.', '15', '+'],
@@ -398,6 +398,22 @@ def test_summary(self):
         self.assertEqual(subprocess.call(command_basic), 0)
         self.assertEqual(subprocess.call(command_full), 0)
 
+    def test_bed2bedgraph(self):
+        command_basic = [
+            'iCount', 'bedgraph', self.cross_links, self.tmp1,
+            '-S', '40',  # Supress lower than ERROR messages.
+        ]
+
+        command_full = [
+            'iCount', 'bedgraph', self.cross_links, self.tmp1,
+            '--name', 'Name.',
+            '--description', 'Description.',
+            '-S', '40',  # Supress lower than ERROR messages.
+        ]
+
+        self.assertEqual(subprocess.call(command_basic), 0)
+        self.assertEqual(subprocess.call(command_full), 0)
+
 
 if __name__ == '__main__':
     unittest.main()

iCount/tests/test_files.py

@@ -178,5 +178,50 @@ def test_fastq_file_write(self):
         self.assertEqual(result, expected)
 
 
+class TestBedGraph(unittest.TestCase):
+
+    def setUp(self):
+        warnings.simplefilter("ignore", ResourceWarning)
+
+        bed_data = [
+            ['1', '4', '5', '.', '5', '+'],
+            ['1', '5', '6', '.', '1', '+'],
+            ['1', '5', '6', '.', '1', '-'],
+            ['2', '5', '6', '.', '3', '+'],
+        ]
+        self.bed = make_file_from_list(bed_data, extension='bed')
+        self.bedgraph = get_temp_file_name(extension='bedgraph')
+
+    def test_bed2bedgraph(self):
+        iCount.files.bedgraph.bed2bedgraph(self.bed, self.bedgraph)
+        expected = [
+            ['track type=bedGraph name="User Track" description="User Supplied Track"'],
+            ['1', '4', '5', '+5'],
+            ['1', '5', '6', '+1'],
+            ['1', '5', '6', '-1'],
+            ['2', '5', '6', '+3'],
+        ]
+        result = make_list_from_file(self.bedgraph, fields_separator='\t')
+        self.assertEqual(result, expected)
+
+    def test_bed2bedgraph_params(self):
+        """
+        Test with custom ``name`` and ``description`` parameters.
+
+        Note that ``name`` is too long and is trimmed to 15 characters.
+        """
+        iCount.files.bedgraph.bed2bedgraph(
+            self.bed, self.bedgraph, name='Longer than 15 chars.', description='Custom text.')
+        expected = [
+            ['track type=bedGraph name="Longer than 15 " description="Custom text."'],
+            ['1', '4', '5', '+5'],
+            ['1', '5', '6', '+1'],
+            ['1', '5', '6', '-1'],
+            ['2', '5', '6', '+3'],
+        ]
+        result = make_list_from_file(self.bedgraph, fields_separator='\t')
+        self.assertEqual(result, expected)
+
+
 if __name__ == '__main__':
     unittest.main()

iCount/tests/utils.py

@@ -155,7 +155,7 @@ def make_aligned_segment(data, rnd_seed=None):
     segment.next_reference_start = 0
     segment.template_length = 0
 
-    length = sum([n2 for (n1, n2) in segment.cigar if n1 in [0, 2, 3, 7, 8]])
+    length = sum([n2 for (n1, n2) in segment.cigar if n1 in [0, 1, 4, 7, 8]])
     segment.query_sequence = make_sequence(size=length, include_n=True, rnd_seed=rnd_seed)
     segment.query_qualities = pysam.qualitystring_to_array(
         make_quality_scores(size=length, rnd_seed=rnd_seed))

setup.py

@@ -56,7 +56,7 @@
             'examples/*.sh',
         ]
     },
-    install_requires={
+    install_requires=[
         'numpy',
         'pandas',
         'cutadapt>=1.10',
@@ -65,7 +65,7 @@
         'numpydoc',
         'sphinx>=1.4',
         'matplotlib',
-    },
+    ],
     extras_require={
         'docs': [
             'docutils',