Skip to content

check_integrity

ensembl.io.genomio.manifest.check_integrity

Compare the genomic data in a DNA FASTA file, seq_region JSON, gene models GFF3 and peptide FASTA to ensure their contents are in sync.

IntegrityTool

Check the integrity of sequence and annotation files in the genome

Source code in src/python/ensembl/io/genomio/manifest/check_integrity.py
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
class IntegrityTool:
    """Check the integrity of sequence and annotation files in the genome"""

    def __init__(
        self,
        manifest_file: Path,
        ignore_final_stops: bool = False,
        no_fail: bool = False,
    ) -> None:
        self.manifest = ManifestStats(manifest_file)
        self.ignore_final_stops = False
        self.set_ignore_final_stops(ignore_final_stops)
        self.errors: list[str] = []
        self.no_fail = no_fail

    def add_errors(self, errors: list[str] | str) -> None:
        """Store the given errors (list or single string) in the list of all errors."""
        if isinstance(errors, str):
            self.errors.append(errors)
        else:
            self.errors += errors

    def check_integrity(self) -> None:
        """Load files listed in the manifest.json and check the integrity.
        Check if the files are correct by verifying the MD5 hash.
        Check if translation, functional annotation and sequence region ids
        and lengths are consistent with the information in gff.
        Compare sequence length from fasta_dna file to seq_region.json metadata.
        """

        # Load the manifest integrity counts
        manifest = self.manifest
        manifest.prepare_integrity_data()

        genome = manifest.genome

        dna = manifest.get_lengths("dna_sequences")
        pep = manifest.get_lengths("peptide_sequences")
        seq_lengths = manifest.get_lengths("seq_regions")
        seq_circular = manifest.get_circular("seq_regions")

        agp_seqr = manifest.get_lengths("agp")

        # Then, run the checks
        self._check_genome(genome)

        if self.manifest.errors:
            errors_str = "\n".join(self.manifest.errors)
            raise InvalidIntegrityError(f"Manifest files parsing failed:\n{errors_str}")

        # Check gff3
        if manifest.has_lengths("gff3_genes"):
            gff_genes = manifest.get_lengths("gff3_genes")
            gff_seq_regions = manifest.get_lengths("gff3_seq_regions")
            gff_translations = manifest.get_lengths("gff3_translations")
            gff_all_translations = manifest.get_lengths("gff3_all_translations")
            gff_transposable_elements = manifest.get_lengths("gff3_transposable_elements")

            ann_genes = manifest.get_lengths("ann_genes")
            ann_translations = manifest.get_lengths("ann_translations")
            ann_transposable_elements = manifest.get_lengths("ann_transposable_elements")

            # Check fasta_pep.fa integrity
            # The sequence length and id retrieved from the fasta_pep file
            # and compared to the translated CDS id and length in the gff
            # We do not compare the peptide lengths because of sequence edits
            if pep:
                tr_errors = self.check_lengths(
                    pep, gff_translations, "Fasta translations vs gff", special_diff=True
                )
                if len(tr_errors) > 0:
                    # The pseudo CDSs are included in this check
                    # Pseudo CDSs are not translated, if the pseudo translation ids are not ignored
                    # in the gff it will give an error
                    tr_errors_all = self.check_lengths(
                        pep,
                        gff_all_translations,
                        "Fasta translations vs gff (include pseudo CDS)",
                        special_diff=True,
                    )
                    if tr_errors_all:
                        self.add_errors(tr_errors)
                        self.add_errors(tr_errors_all)

            # Check functional_annotation.json integrity
            # Gene ids, translated CDS ids and translated CDSs
            # including pseudogenes are compared to the gff
            if ann_genes:
                self.add_errors(self.check_ids(ann_genes, gff_genes, "Gene ids metadata vs gff"))
                tr_id_errors = self.check_ids(
                    ann_translations, gff_translations, "Translation ids metadata vs gff"
                )
                if tr_id_errors:
                    tr_id_errors_all = self.check_ids(
                        ann_translations,
                        gff_all_translations,
                        "Translation ids metadata vs gff (include pseudo CDS)",
                    )
                    if tr_id_errors_all:
                        self.add_errors(tr_id_errors)
                        self.add_errors(tr_id_errors_all)
                self.add_errors(
                    self.check_ids(
                        ann_transposable_elements,
                        gff_transposable_elements,
                        "TE ids metadata vs gff",
                    )
                )

            self.check_seq_region_lengths(
                seq_lengths, gff_seq_regions, "seq_regions JSON vs GFF3 lengths", seq_circular
            )

        self.check_seq_region_lengths(seq_lengths, dna, "seq_regions JSON vs DNA lengths")
        self.check_seq_region_lengths(seq_lengths, agp_seqr, "seq_regions JSON vs AGPs lengths")

        if self.errors:
            errors_str = "\n".join(self.errors)
            message = f"Integrity test failed:\n{errors_str}"
            if self.no_fail:
                print(message)
            else:
                raise InvalidIntegrityError(message)

    def set_ignore_final_stops(self, ignore_final_stops: bool) -> None:
        """Set ignore_final_stops (when calculating peptide length) for this tool and the manifest."""
        self.ignore_final_stops = ignore_final_stops
        self.manifest.ignore_final_stops = ignore_final_stops

    def _check_genome(self, genome: dict[str, Any]) -> None:
        """Check if the accession is correct in genome.json."""
        genome_accession = genome.get("assembly", {}).get("accession", "")
        if not genome_accession:
            return
        if not re.match(r"GC[AF]_\d{9}(\.\d+)?", genome_accession):
            self.add_errors(f"Genome assembly accession is wrong: '{genome_accession}'")

    def check_ids(self, list1: dict[str, Any], list2: dict[str, Any], name: str) -> list[str]:
        """Compare the ids in list1 and list2.

        Args:
            list1: Sequence IDs retrieved from `functional.json`.
            list2: Sequence IDs retrieved from the GFF3 file.
            name: Source name.

        Return:
            List of message errors of sequence IDs found only in one of the lists provided.
        """

        only1 = []
        only2 = []
        common = []

        for item_id in list1:
            if item_id in list2:
                common.append(item_id)
            else:
                only1.append(item_id)
        for item_id in list2:
            if item_id not in common:
                only2.append(item_id)

        errors = []
        if common:
            logging.info(f"{len(common)} common elements in {name}")
        if only1:
            errors.append(f"{len(only1)} only in first list in {name} (first: {only1[0]})")
            logging.debug(f"{len(only1)} only in first list in {name}")
        if only2:
            errors.append(f"{len(only2)} only in second list in {name} (first: {only2[0]})")
            logging.debug(f"{len(only1)} only in second list in {name}")

        return errors

    def check_lengths(
        self,
        list1: dict[str, int],
        list2: dict[str, int],
        name: str,
        *,
        allowed_len_diff: int | None = None,
        special_diff: bool = False,
    ) -> list[str]:
        """Check the difference in ids and length between list1 and list2.
            There are a few special cases here where we allow a certain asymmetry
            by changing the values of the arguments.

        Args:
            list1: dict containing length and id of the sequence from fasta files.
            list2: dict containing length and id in the retrieved from the gff.
            name:  string

        allowed_len_diff : None to to not accept differences in length between list1 and list2.
            The value can be changed based on how much difference in sequence length we are wanting to accept.

        special_diff: set as False when no special length difference is expected between the lists.
                    This can be changed if we want to report common sequences with 1 BP difference.

        Returns:
            Error if there is a difference in length or ids between the lists.
        """

        # check list differences, checks if abs(values diff) < allowed_len_diff

        set1 = frozenset(list1)
        set2 = frozenset(list2)
        list1_2 = list(set1 - set2)
        list2_1 = list(set2 - set1)

        errors = []
        if len(list1_2) > 0:
            errors.append(f"{name}: {len(list1_2)} from the first list only (i.e. {list1_2[0]})")
        if len(list2_1) > 0:
            errors.append(f"{name}: {len(list2_1)} from the second list only (i.e. {list2_1[0]})")

        common_len = 0
        if allowed_len_diff is None:
            common_len = len(set1 & set2)
        else:
            # check for the sequence length difference
            diff_len_list: list[str] = []
            diff_len_special_list: list[str] = []
            for e in set1 & set2:
                dl12 = list1[e] - list2[e]
                if abs(dl12) <= allowed_len_diff:
                    common_len += 1
                else:
                    _dlist = diff_len_list
                    # Special case: 1 AA /BP shorter,
                    #   so assuming the stop codon is not included in the CDS (when it should be)
                    if dl12 == 1 and special_diff:
                        _dlist = diff_len_special_list
                    _dlist.append(f"{e}: {list1[e]}, {list2[e]}")
            if diff_len_special_list:
                errors.append(
                    (
                        f"{len(diff_len_special_list)} common elements with one BP/AA length diff for {name}"
                        f"(e.g. {diff_len_special_list[0]})"
                    )
                )
            if diff_len_list:
                errors.append(
                    (
                        f"{len(diff_len_list)} common elements with length diff for {name}"
                        f"(e.g. {diff_len_list[0]})"
                    )
                )
        if common_len > 0:
            logging.warning(f"{common_len} common elements between lists for {name}")

        return errors

    def check_seq_region_lengths(
        self,
        seqrs: dict[str, Any] | None,
        feats: dict[str, Any] | None,
        name: str,
        circular: dict[str, Any] | None = None,
    ) -> None:
        """Check the integrity of seq_region.json file by comparing the length of the sequence
            to fasta files and the gff.

            Seq_region file is in json format containing the metadata of the sequence.
            It contains sequence id, length, location and the synonyms for the sequence name
            from different sources.

        Args:
            seqs: Sequence name and length retrieved from seq_region.json file.
            feats: Sequence name and length retrieved from the fasta and gff file.
            name: Name of the check to show in the logs.
            circular: Whether any sequence is circular.

        Returns:
            Error if there are common sequences with difference in ids
            and if the sequences are not consistent in the files.
        """
        if not seqrs or not feats:
            return
        comp = self._compare_seqs(seqrs, feats, circular)

        common = comp["common"]
        diff = comp["diff"]
        diff_circular = comp["diff_circular"]
        only_seqr = comp["only_seqr"]
        only_feat = comp["only_feat"]

        if common:
            logging.info(f"{len(common)} common elements in {name}")
        if diff_circular:
            example = diff_circular[0]
            logging.info(f"{len(diff_circular)} differences for circular elements in {name} (e.g. {example})")
        if diff:
            self.add_errors(f"{len(diff)} common elements with higher length in {name} (e.g. {diff[0]})")
        if only_seqr:
            # Not an error!
            logging.info(f"{len(only_seqr)} only in seq_region list in {name} (first: {only_seqr[0]})")
        if only_feat:
            self.add_errors(f"{len(only_feat)} only in second list in {name} (first: {only_feat[0]})")

    def _compare_seqs(
        self, seqrs: dict[str, Any], feats: dict[str, Any], circular: dict[str, Any] | None = None
    ) -> dict[str, list[str]]:
        """Give the intersection and other comparison between two groups of sequences.

        Args:
            seqs: Sequence name and length retrieved from seq_region.json file.
            feats: Sequence name and length retrieved from the fasta and gff file.
            circular: Whether any sequence is circular.

        Returns: Dict with 5 stats:
            common: Common elements.
            only_seqr: Elements only in the first one.
            only_feat: Elements only in the second one.
            diff: Elements that differ.
            diff_circular: Elements that differ in a circular sequence.

        """
        comp: dict[str, list[str]] = {
            "common": [],
            "only_seqr": [],
            "only_feat": [],
            "diff": [],
            "diff_circular": [],
        }

        for seq_id in seqrs:
            if seq_id in feats:
                # Check that feature is within the seq_region length
                if feats[seq_id] > seqrs[seq_id]:
                    diff_str = f"{seq_id}: {seqrs[seq_id]} vs {feats[seq_id]}"
                    if circular and circular.get(seq_id, False):
                        comp["diff_circular"].append(diff_str)
                    else:
                        comp["diff"].append(diff_str)
                else:
                    comp["common"].append(seq_id)
            else:
                comp["only_seqr"].append(seq_id)

        for seq_id in feats:
            if (
                seq_id not in comp["common"]
                and seq_id not in comp["diff"]
                and seq_id not in comp["diff_circular"]
                and seq_id not in seqrs
            ):
                comp["only_feat"].append(seq_id)

        return comp

errors = [] instance-attribute

ignore_final_stops = False instance-attribute

manifest = ManifestStats(manifest_file) instance-attribute

no_fail = no_fail instance-attribute

add_errors(errors)

Store the given errors (list or single string) in the list of all errors.

Source code in src/python/ensembl/io/genomio/manifest/check_integrity.py
47
48
49
50
51
52
def add_errors(self, errors: list[str] | str) -> None:
    """Store the given errors (list or single string) in the list of all errors."""
    if isinstance(errors, str):
        self.errors.append(errors)
    else:
        self.errors += errors

check_ids(list1, list2, name)

Compare the ids in list1 and list2.

Parameters:

Name Type Description Default
list1 dict[str, Any]

Sequence IDs retrieved from functional.json.

required
list2 dict[str, Any]

Sequence IDs retrieved from the GFF3 file.

required
name str

Source name.

required
Return

List of message errors of sequence IDs found only in one of the lists provided.

Source code in src/python/ensembl/io/genomio/manifest/check_integrity.py
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
def check_ids(self, list1: dict[str, Any], list2: dict[str, Any], name: str) -> list[str]:
    """Compare the ids in list1 and list2.

    Args:
        list1: Sequence IDs retrieved from `functional.json`.
        list2: Sequence IDs retrieved from the GFF3 file.
        name: Source name.

    Return:
        List of message errors of sequence IDs found only in one of the lists provided.
    """

    only1 = []
    only2 = []
    common = []

    for item_id in list1:
        if item_id in list2:
            common.append(item_id)
        else:
            only1.append(item_id)
    for item_id in list2:
        if item_id not in common:
            only2.append(item_id)

    errors = []
    if common:
        logging.info(f"{len(common)} common elements in {name}")
    if only1:
        errors.append(f"{len(only1)} only in first list in {name} (first: {only1[0]})")
        logging.debug(f"{len(only1)} only in first list in {name}")
    if only2:
        errors.append(f"{len(only2)} only in second list in {name} (first: {only2[0]})")
        logging.debug(f"{len(only1)} only in second list in {name}")

    return errors

check_integrity()

Load files listed in the manifest.json and check the integrity. Check if the files are correct by verifying the MD5 hash. Check if translation, functional annotation and sequence region ids and lengths are consistent with the information in gff. Compare sequence length from fasta_dna file to seq_region.json metadata.

Source code in src/python/ensembl/io/genomio/manifest/check_integrity.py
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
def check_integrity(self) -> None:
    """Load files listed in the manifest.json and check the integrity.
    Check if the files are correct by verifying the MD5 hash.
    Check if translation, functional annotation and sequence region ids
    and lengths are consistent with the information in gff.
    Compare sequence length from fasta_dna file to seq_region.json metadata.
    """

    # Load the manifest integrity counts
    manifest = self.manifest
    manifest.prepare_integrity_data()

    genome = manifest.genome

    dna = manifest.get_lengths("dna_sequences")
    pep = manifest.get_lengths("peptide_sequences")
    seq_lengths = manifest.get_lengths("seq_regions")
    seq_circular = manifest.get_circular("seq_regions")

    agp_seqr = manifest.get_lengths("agp")

    # Then, run the checks
    self._check_genome(genome)

    if self.manifest.errors:
        errors_str = "\n".join(self.manifest.errors)
        raise InvalidIntegrityError(f"Manifest files parsing failed:\n{errors_str}")

    # Check gff3
    if manifest.has_lengths("gff3_genes"):
        gff_genes = manifest.get_lengths("gff3_genes")
        gff_seq_regions = manifest.get_lengths("gff3_seq_regions")
        gff_translations = manifest.get_lengths("gff3_translations")
        gff_all_translations = manifest.get_lengths("gff3_all_translations")
        gff_transposable_elements = manifest.get_lengths("gff3_transposable_elements")

        ann_genes = manifest.get_lengths("ann_genes")
        ann_translations = manifest.get_lengths("ann_translations")
        ann_transposable_elements = manifest.get_lengths("ann_transposable_elements")

        # Check fasta_pep.fa integrity
        # The sequence length and id retrieved from the fasta_pep file
        # and compared to the translated CDS id and length in the gff
        # We do not compare the peptide lengths because of sequence edits
        if pep:
            tr_errors = self.check_lengths(
                pep, gff_translations, "Fasta translations vs gff", special_diff=True
            )
            if len(tr_errors) > 0:
                # The pseudo CDSs are included in this check
                # Pseudo CDSs are not translated, if the pseudo translation ids are not ignored
                # in the gff it will give an error
                tr_errors_all = self.check_lengths(
                    pep,
                    gff_all_translations,
                    "Fasta translations vs gff (include pseudo CDS)",
                    special_diff=True,
                )
                if tr_errors_all:
                    self.add_errors(tr_errors)
                    self.add_errors(tr_errors_all)

        # Check functional_annotation.json integrity
        # Gene ids, translated CDS ids and translated CDSs
        # including pseudogenes are compared to the gff
        if ann_genes:
            self.add_errors(self.check_ids(ann_genes, gff_genes, "Gene ids metadata vs gff"))
            tr_id_errors = self.check_ids(
                ann_translations, gff_translations, "Translation ids metadata vs gff"
            )
            if tr_id_errors:
                tr_id_errors_all = self.check_ids(
                    ann_translations,
                    gff_all_translations,
                    "Translation ids metadata vs gff (include pseudo CDS)",
                )
                if tr_id_errors_all:
                    self.add_errors(tr_id_errors)
                    self.add_errors(tr_id_errors_all)
            self.add_errors(
                self.check_ids(
                    ann_transposable_elements,
                    gff_transposable_elements,
                    "TE ids metadata vs gff",
                )
            )

        self.check_seq_region_lengths(
            seq_lengths, gff_seq_regions, "seq_regions JSON vs GFF3 lengths", seq_circular
        )

    self.check_seq_region_lengths(seq_lengths, dna, "seq_regions JSON vs DNA lengths")
    self.check_seq_region_lengths(seq_lengths, agp_seqr, "seq_regions JSON vs AGPs lengths")

    if self.errors:
        errors_str = "\n".join(self.errors)
        message = f"Integrity test failed:\n{errors_str}"
        if self.no_fail:
            print(message)
        else:
            raise InvalidIntegrityError(message)

check_lengths(list1, list2, name, *, allowed_len_diff=None, special_diff=False)

Check the difference in ids and length between list1 and list2. There are a few special cases here where we allow a certain asymmetry by changing the values of the arguments.

Parameters:

Name Type Description Default
list1 dict[str, int]

dict containing length and id of the sequence from fasta files.

required
list2 dict[str, int]

dict containing length and id in the retrieved from the gff.

required
name str

string

required
None to to not accept differences in length between list1 and list2.

The value can be changed based on how much difference in sequence length we are wanting to accept.

set as False when no special length difference is expected between the lists.

This can be changed if we want to report common sequences with 1 BP difference.

Returns:

Type Description
list[str]

Error if there is a difference in length or ids between the lists.

Source code in src/python/ensembl/io/genomio/manifest/check_integrity.py
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
def check_lengths(
    self,
    list1: dict[str, int],
    list2: dict[str, int],
    name: str,
    *,
    allowed_len_diff: int | None = None,
    special_diff: bool = False,
) -> list[str]:
    """Check the difference in ids and length between list1 and list2.
        There are a few special cases here where we allow a certain asymmetry
        by changing the values of the arguments.

    Args:
        list1: dict containing length and id of the sequence from fasta files.
        list2: dict containing length and id in the retrieved from the gff.
        name:  string

    allowed_len_diff : None to to not accept differences in length between list1 and list2.
        The value can be changed based on how much difference in sequence length we are wanting to accept.

    special_diff: set as False when no special length difference is expected between the lists.
                This can be changed if we want to report common sequences with 1 BP difference.

    Returns:
        Error if there is a difference in length or ids between the lists.
    """

    # check list differences, checks if abs(values diff) < allowed_len_diff

    set1 = frozenset(list1)
    set2 = frozenset(list2)
    list1_2 = list(set1 - set2)
    list2_1 = list(set2 - set1)

    errors = []
    if len(list1_2) > 0:
        errors.append(f"{name}: {len(list1_2)} from the first list only (i.e. {list1_2[0]})")
    if len(list2_1) > 0:
        errors.append(f"{name}: {len(list2_1)} from the second list only (i.e. {list2_1[0]})")

    common_len = 0
    if allowed_len_diff is None:
        common_len = len(set1 & set2)
    else:
        # check for the sequence length difference
        diff_len_list: list[str] = []
        diff_len_special_list: list[str] = []
        for e in set1 & set2:
            dl12 = list1[e] - list2[e]
            if abs(dl12) <= allowed_len_diff:
                common_len += 1
            else:
                _dlist = diff_len_list
                # Special case: 1 AA /BP shorter,
                #   so assuming the stop codon is not included in the CDS (when it should be)
                if dl12 == 1 and special_diff:
                    _dlist = diff_len_special_list
                _dlist.append(f"{e}: {list1[e]}, {list2[e]}")
        if diff_len_special_list:
            errors.append(
                (
                    f"{len(diff_len_special_list)} common elements with one BP/AA length diff for {name}"
                    f"(e.g. {diff_len_special_list[0]})"
                )
            )
        if diff_len_list:
            errors.append(
                (
                    f"{len(diff_len_list)} common elements with length diff for {name}"
                    f"(e.g. {diff_len_list[0]})"
                )
            )
    if common_len > 0:
        logging.warning(f"{common_len} common elements between lists for {name}")

    return errors

check_seq_region_lengths(seqrs, feats, name, circular=None)

Check the integrity of seq_region.json file by comparing the length of the sequence to fasta files and the gff.

Seq_region file is in json format containing the metadata of the sequence.
It contains sequence id, length, location and the synonyms for the sequence name
from different sources.

Parameters:

Name Type Description Default
seqs

Sequence name and length retrieved from seq_region.json file.

required
feats dict[str, Any] | None

Sequence name and length retrieved from the fasta and gff file.

required
name str

Name of the check to show in the logs.

required
circular dict[str, Any] | None

Whether any sequence is circular.

None

Returns:

Type Description
None

Error if there are common sequences with difference in ids

None

and if the sequences are not consistent in the files.

Source code in src/python/ensembl/io/genomio/manifest/check_integrity.py
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
def check_seq_region_lengths(
    self,
    seqrs: dict[str, Any] | None,
    feats: dict[str, Any] | None,
    name: str,
    circular: dict[str, Any] | None = None,
) -> None:
    """Check the integrity of seq_region.json file by comparing the length of the sequence
        to fasta files and the gff.

        Seq_region file is in json format containing the metadata of the sequence.
        It contains sequence id, length, location and the synonyms for the sequence name
        from different sources.

    Args:
        seqs: Sequence name and length retrieved from seq_region.json file.
        feats: Sequence name and length retrieved from the fasta and gff file.
        name: Name of the check to show in the logs.
        circular: Whether any sequence is circular.

    Returns:
        Error if there are common sequences with difference in ids
        and if the sequences are not consistent in the files.
    """
    if not seqrs or not feats:
        return
    comp = self._compare_seqs(seqrs, feats, circular)

    common = comp["common"]
    diff = comp["diff"]
    diff_circular = comp["diff_circular"]
    only_seqr = comp["only_seqr"]
    only_feat = comp["only_feat"]

    if common:
        logging.info(f"{len(common)} common elements in {name}")
    if diff_circular:
        example = diff_circular[0]
        logging.info(f"{len(diff_circular)} differences for circular elements in {name} (e.g. {example})")
    if diff:
        self.add_errors(f"{len(diff)} common elements with higher length in {name} (e.g. {diff[0]})")
    if only_seqr:
        # Not an error!
        logging.info(f"{len(only_seqr)} only in seq_region list in {name} (first: {only_seqr[0]})")
    if only_feat:
        self.add_errors(f"{len(only_feat)} only in second list in {name} (first: {only_feat[0]})")

set_ignore_final_stops(ignore_final_stops)

Set ignore_final_stops (when calculating peptide length) for this tool and the manifest.

Source code in src/python/ensembl/io/genomio/manifest/check_integrity.py
156
157
158
159
def set_ignore_final_stops(self, ignore_final_stops: bool) -> None:
    """Set ignore_final_stops (when calculating peptide length) for this tool and the manifest."""
    self.ignore_final_stops = ignore_final_stops
    self.manifest.ignore_final_stops = ignore_final_stops

main()

Main entrypoint.

Source code in src/python/ensembl/io/genomio/manifest/check_integrity.py
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
def main() -> None:
    """Main entrypoint."""
    parser = ArgumentParser(description=__doc__)
    parser.add_argument_src_path("--manifest_file", required=True, help="Manifest file for the data to check")
    parser.add_argument(
        "--ignore_final_stops", action="store_true", help="Ignore final stop when calculating peptide length"
    )
    parser.add_argument(
        "--no_fail", action="store_true", help="In case of errors, don't fail but print errors to stdout."
    )
    parser.add_argument("--version", action="version", version=ensembl.io.genomio.__version__)
    parser.add_log_arguments(add_log_file=True)
    args = parser.parse_args()
    init_logging_with_args(args)

    inspector = IntegrityTool(args.manifest_file, args.ignore_final_stops, args.no_fail)
    inspector.check_integrity()