feat (moPepGen): Added tally table to parseREDITools, parseCIRCExplorer, and parseRMATS

Author: zhuchcn

CHANGELOG.md

@@ -14,6 +14,8 @@ This project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.htm
 
 - Added --skip-failed flag to callVariant, parseArriba, parserSTARFusion, parseFusionCatcher.
 
+- Added tally table to parseREDITools, parseCIRCExplorer, and parseRMATS
+
 ## [1.4.6-rc2] - 2025-03-03
 
 ### Fixed

moPepGen/cli/parse_circexplorer.py

@@ -4,14 +4,18 @@
 [callVariant](call-variant.md). Noted that only known circRNA is supported (
 \*_circular_known.txt) """
 from __future__ import annotations
+from typing import TYPE_CHECKING
 import argparse
-from typing import List, Dict
 from pathlib import Path
 from moPepGen import get_logger, circ, err
 from moPepGen.parser import CIRCexplorerParser
 from moPepGen.cli import common
 
 
+if TYPE_CHECKING:
+    from typing import List, Dict
+    from logging import Logger
+
 INPUT_FILE_FORMATS = ['.tsv', '.txt']
 OUTPUT_FILE_FORMATS = ['.gvf']
 
@@ -74,16 +78,48 @@ def add_subparser_parse_circexplorer(subparsers:argparse._SubParsersAction):
         default='-100,5',
         metavar='<number>'
     )
+    common.add_args_skip_failed(p)
     common.add_args_source(p)
     common.add_args_reference(p, genome=False, proteome=False)
     common.add_args_debug_level(p)
     p.set_defaults(func=parse_circexplorer)
     common.print_help_if_missing_args(p)
     return p
 
+class TallyTable():
+    """ Tally table """
+    def __init__(self, logger:Logger):
+        """ Constructor """
+        self.total:int = 0
+        self.succeed:int = 0
+        self.skipped:TallyTableSkipped = TallyTableSkipped()
+        self.logger = logger
+
+    def log(self):
+        """ Show tally results """
+        self.logger.info("Summary:")
+        self.logger.info("Totally records read: %i", self.total)
+        self.logger.info("Records successfully processed: %i", self.succeed)
+        self.logger.info("Records skipped: %i", self.skipped.total)
+        if self.skipped.total > 0:
+            self.logger.info("Out of those skipped,")
+            self.logger.info("    Invalid circRNA record: %i", self.skipped.invalid_record)
+            self.logger.info("    Insufficient evidence: %i", self.skipped.insufficient_evidence)
+
+class TallyTableSkipped():
+    """ Tally table for failed ones """
+    def __init__(self):
+        """ constructor """
+        self.invalid_gene_id:int = 0
+        self.invalid_position:int = 0
+        self.insufficient_evidence:int = 0
+        self.invalid_record:int = 0
+        self.total:int = 0
+
 def parse_circexplorer(args:argparse.Namespace):
     """ Parse circexplorer known circRNA results. """
     logger = get_logger()
+    tally = TallyTable(logger)
 
     input_path:Path = args.input_path
     output_path:Path = args.output_path
@@ -104,11 +140,16 @@ def parse_circexplorer(args:argparse.Namespace):
     circ_records:Dict[str, List[circ.CircRNAModel]] = {}
 
     for record in CIRCexplorerParser.parse(input_path, args.circexplorer3):
+        tally.total += 1
         if not args.circexplorer3:
             if not record.is_valid(args.min_read_number):
+                tally.skipped.total += 1
+                tally.skipped.insufficient_evidence += 1
                 continue
         elif not record.is_valid(args.min_read_number, args.min_fbr_circ, \
                 args.min_circ_score):
+            tally.skipped.total += 1
+            tally.skipped.insufficient_evidence += 1
             continue
         try:
             circ_record = record.convert_to_circ_rna(anno, intron_start_range,
@@ -119,13 +160,17 @@ def parse_circexplorer(args:argparse.Namespace):
                 " Skipping it from parsing.",
                 record.name, record.isoform_name
             )
+            tally.skipped.invalid_record += 1
+            tally.skipped.total += 1
             continue
         except err.IntronNotFoundError:
             logger.warning(
                 "The CIRCexplorer record %s from transcript %s contains an unknown"
                 " intron. Skipping it from parsing.",
                 record.name, record.isoform_name
             )
+            tally.skipped.invalid_record += 1
+            tally.skipped.total += 1
             continue
         except:
             logger.error('Exception raised from record: %s', record.name)
@@ -135,21 +180,20 @@ def parse_circexplorer(args:argparse.Namespace):
             circ_records[gene_id] = []
         circ_records[gene_id].append(circ_record)
 
-    if not circ_records:
-        logger.warning('No variant record is saved.')
-        return
+    if circ_records:
+        genes_rank = anno.get_genes_rank()
+        ordered_keys = sorted(circ_records.keys(), key=lambda x:genes_rank[x])
 
-    genes_rank = anno.get_genes_rank()
-    ordered_keys = sorted(circ_records.keys(), key=lambda x:genes_rank[x])
+        records = []
+        for key in ordered_keys:
+            val = circ_records[key]
+            records.extend(val)
 
-    records = []
-    for key in ordered_keys:
-        val = circ_records[key]
-        records.extend(val)
+        metadata = common.generate_metadata(args)
 
-    metadata = common.generate_metadata(args)
+        with open(output_path, 'w') as handle:
+            circ.io.write(records, metadata, handle)
 
-    with open(output_path, 'w') as handle:
-        circ.io.write(records, metadata, handle)
+        logger.info("CircRNA records written to disk.")
 
-    logger.info("CircRNA records written to disk.")
+    tally.log()

moPepGen/cli/parse_reditools.py

@@ -4,13 +4,17 @@
 [callVariant](call-variant.md)
 """
 from __future__ import annotations
+from typing import TYPE_CHECKING
 import argparse
 from pathlib import Path
-from typing import Dict, List
 from moPepGen import get_logger, seqvar, parser
 from moPepGen.cli import common
 
 
+if TYPE_CHECKING:
+    from typing import Dict, List
+    from logging import Logger
+
 INPUT_FILE_FORMATS = ['.tsv', '.txt']
 OUTPUT_FILE_FORMATS = ['.gvf']
 
@@ -76,9 +80,26 @@ def add_subparser_parse_reditools(subparsers:argparse._SubParsersAction):
     common.print_help_if_missing_args(p)
     return p
 
+class TallyTable():
+    """ Tally table """
+    def __init__(self, logger:Logger):
+        """ Constructor """
+        self.total:int = 0
+        self.succeed:int = 0
+        self.skipped:int = 0
+        self.logger = logger
+
+    def log(self):
+        """ Show tally results """
+        self.logger.info("Summary:")
+        self.logger.info("Totally records read: %i", self.total)
+        self.logger.info("Records successfully processed: %i", self.succeed)
+        self.logger.info("Records skipped: %i", self.skipped)
+
 def parse_reditools(args:argparse.Namespace) -> None:
     """ Parse REDItools output and save it in the GVF format. """
     logger = get_logger()
+    tally = TallyTable(logger)
     # unpack args
     table_file:Path = args.input_path
     output_path:Path = args.output_path
@@ -102,13 +123,18 @@ def parse_reditools(args:argparse.Namespace) -> None:
     variants:Dict[str, List[seqvar.VariantRecord]] = {}
 
     for record in parser.REDItoolsParser.parse(table_file, transcript_id_column):
+        tally.total += 1
         _vars = record.convert_to_variant_records(
             anno=anno,
             min_coverage_alt=min_coverage_alt,
             min_frequency_alt=min_frequency_alt,
             min_coverage_rna=min_coverage_rna,
             min_coverage_dna=min_coverage_dna
         )
+        if not _vars:
+            tally.skipped += 1
+        else:
+            tally.succeed += 1
         for variant in _vars:
             gene_id = variant.location.seqname
             if gene_id not in variants:
@@ -139,3 +165,5 @@ def parse_reditools(args:argparse.Namespace) -> None:
     seqvar.io.write(all_records, output_path, metadata)
 
     logger.info('Variants written to disk.')
+
+    tally.log()

moPepGen/cli/parse_rmats.py

@@ -7,14 +7,18 @@
 [callVariant](call-variant.md)
 """
 from __future__ import annotations
+from typing import TYPE_CHECKING
 import argparse
-from typing import Dict, Set
 from pathlib import Path
 from moPepGen import get_logger, seqvar
 from moPepGen.parser import RMATSParser
 from moPepGen.cli import common
 
 
+if TYPE_CHECKING:
+    from typing import Dict, Set
+    from logging import Logger
+
 INPUT_FILE_FORMATS = ['.tsv', '.txt']
 OUTPUT_FILE_FORMATS = ['.gvf']
 
@@ -102,9 +106,27 @@ def add_subparser_parse_rmats(subparsers:argparse._SubParsersAction):
     common.print_help_if_missing_args(p)
     return p
 
+class TallyTable():
+    """ Tally table """
+    def __init__(self, logger:Logger):
+        """ Constructor """
+        self.total:int = 0
+        self.succeed:int = 0
+        self.skipped:int = 0
+        self.logger = logger
+
+    def log(self):
+        """ Show tally results """
+        self.logger.info("Summary:")
+        self.logger.info("Totally records read: %i", self.total)
+        self.logger.info("Records successfully processed: %i", self.succeed)
+        self.logger.info("Records skipped: %i", self.skipped)
+
+
 def parse_rmats(args:argparse.Namespace) -> None:
     """ Parse rMATS results into TSV """
     logger = get_logger()
+    tally = TallyTable(logger)
 
     skipped_exon = args.skipped_exon
     alternative_5 = args.alternative_5_splicing
@@ -136,6 +158,7 @@ def parse_rmats(args:argparse.Namespace) -> None:
         if path:
             logger.info("Start parsing %s file %s", event_type, path)
             for record in RMATSParser.parse(path, event_type):
+                tally.total += 1
                 try:
                     var_records = record.convert_to_variant_records(
                         anno=anno, genome=genome,
@@ -144,6 +167,10 @@ def parse_rmats(args:argparse.Namespace) -> None:
                 except:
                     logger.error(record.gene_id)
                     raise
+                if var_records:
+                    tally.succeed += 1
+                else:
+                    tally.skipped += 1
                 for var_record in var_records:
                     tx_id = var_record.transcript_id
                     if tx_id not in variants:
@@ -168,3 +195,5 @@ def parse_rmats(args:argparse.Namespace) -> None:
     seqvar.io.write(variants_sorted, output_path, metadata)
 
     logger.info('Variants written to disk.')
+
+    tally.log()

test/integration/test_parse_arriba.py

@@ -2,7 +2,7 @@
 import argparse
 import subprocess as sp
 import sys
-from unittest.mock import Mock
+from unittest import mock
 from test.integration import TestCaseIntegration
 from moPepGen import cli
 
@@ -55,12 +55,12 @@ def test_parse_arriba_cli(self):
             print(res.stderr.decode('utf-8'))
             raise
 
+    @mock.patch(
+        "moPepGen.parser.ArribaParser.ArribaRecord.convert_to_variant_records",
+        new=mock.MagicMock(side_effect=ValueError())
+    )
     def test_parse_arriba_skip_failed(self):
         """ Test parseArriba with skip failed """
-        from moPepGen import parser
-        parser.ArribaParser.ArribaRecord.convert_to_variant_records = Mock(
-            side_effect=ValueError()
-        )
         args = self.create_base_args()
         with self.assertRaises(ValueError):
             cli.parse_arriba(args)

test/integration/test_parse_fusion_catcher.py

@@ -3,7 +3,7 @@
 from pathlib import Path
 import subprocess as sp
 import sys
-from unittest.mock import Mock
+from unittest import mock
 from test.unit import load_references
 from test.integration import TestCaseIntegration
 from moPepGen import cli, parser
@@ -81,13 +81,12 @@ def test_parse_fusion_catcher(self):
         self.assertEqual(files, expected)
         self.assert_gvf_order(args.output_path, args.annotation_gtf)
 
+    @mock.patch(
+        "moPepGen.parser.FusionCatcherParser.FusionCatcherRecord.convert_to_variant_records",
+        new=mock.MagicMock(side_effect=ValueError())
+    )
     def test_parse_fusion_catcher_skip_failed(self):
         """ Test parseFusionCatcher with --skip-failed """
-        from moPepGen import parser
-        parser.FusionCatcherParser.FusionCatcherRecord.convert_to_variant_records = Mock(
-            side_effect=ValueError()
-        )
-
         args = self.create_base_args()
         args.input_path = self.data_dir/'fusion/fusion_catcher.txt'
         with self.assertRaises(ValueError):

test/integration/test_parse_star_fusion.py

@@ -2,7 +2,7 @@
 import argparse
 import subprocess as sp
 import sys
-from unittest.mock import Mock
+from unittest import mock
 from test.integration import TestCaseIntegration
 from moPepGen import cli, seqvar
 from moPepGen.cli.common import load_references
@@ -79,12 +79,12 @@ def test_parse_star_fusion_case1(self):
         self.assertEqual(files, expected)
         self.assert_gvf_order(args.output_path, args.annotation_gtf)
 
+    @mock.patch(
+        "moPepGen.parser.STARFusionParser.STARFusionRecord.convert_to_variant_records",
+        new=mock.MagicMock(side_effect=ValueError())
+    )
     def test_parse_star_fusion_skip_failed(self):
         """ test parseSTARFusion case1 """
-        from moPepGen import parser
-        parser.STARFusionParser.STARFusionRecord.convert_to_variant_records = Mock(
-            side_effect=ValueError()
-        )
         args = self.create_base_args()
         args.input_path = self.data_dir/'fusion/star_fusion.txt'
         with self.assertRaises(ValueError):