comparison intersect_and_adjacent.py @ 5:00aab5199327 draft default tip

planemo upload commit 6b68423e1d9bbced95ecfb92e26329a7e8d7df47

4:3e3b5ba626b9

            # Used to map ids back to features later
            rec_a_map = {f.id: f for f in rec_a_i.features}
            rec_b_map = {f.id: f for f in rec_b_i.features}

            rec_a_hits_in_b = []
            rec_b_hits_in_a = []

            for feature in rec_a_i.features:
                if feature.type == "remark" or feature.type == "annotation":
                    continue

                if not stranding:
                    hits = tree_b[
                        int(feature.location.start) : int(feature.location.end)
                    ]
                    for hit in hits:
                        rec_a_hits_in_b.append(rec_b_map[hit.data])
                else:
                    if feature.strand > 0:
                        hits = tree_b_pos[
                            int(feature.location.start) : int(feature.location.end)
                        ]
                    else:
                        hits = tree_b_neg[
                            int(feature.location.start) : int(feature.location.end)
                        ]
                    for hit in hits:
                        rec_a_hits_in_b.append(rec_b_map[hit.data])

            for feature in rec_b_i.features:
                if feature.type == "remark" or feature.type == "annotation":
                    continue

                if not stranding:
                    hits = tree_a[
                        int(feature.location.start) : int(feature.location.end)
                    ]
                    for hit in hits:
                        rec_b_hits_in_a.append(rec_a_map[hit.data])
                else:
                    if feature.strand > 0:
                        hits = tree_a_pos[
                            int(feature.location.start) : int(feature.location.end)
                        ]
                    else:
                        hits = tree_a_neg[
                            int(feature.location.start) : int(feature.location.end)
                        ]
                    for hit in hits:
                        rec_b_hits_in_a.append(rec_a_map[hit.data])

            # Remove duplicate features using sets
            rec_a_out.append(
                SeqRecord(
                    rec_a[iterate].seq,
                    rec_a[iterate].id,
                    rec_a[iterate].name,
                    rec_a[iterate].description,
                    rec_a[iterate].dbxrefs,
                    sorted(set(rec_a_hits_in_b), key=lambda feat: feat.location.start),
                    rec_a[iterate].annotations,
                )
            )
            rec_b_out.append(
                SeqRecord(
                    rec_b[iterate].seq,
                    rec_b[iterate].id,
                    rec_b[iterate].name,
                    rec_b[iterate].description,
                    rec_b[iterate].dbxrefs,
                    sorted(set(rec_b_hits_in_a), key=lambda feat: feat.location.start),
                    rec_b[iterate].annotations,
                )
            )
            iterate += 1

5:00aab5199327

            # Used to map ids back to features later
            rec_a_map = {f.id: f for f in rec_a_i.features}
            rec_b_map = {f.id: f for f in rec_b_i.features}

            rec_a_hits_in_b = {}
            rec_b_hits_in_a = {}

            for feature in rec_a_i.features:
                if feature.type == "remark" or feature.type == "annotation":
                    continue

                if not stranding:
                    hits = tree_b[
                        int(feature.location.start) : int(feature.location.end)
                    ]
                    for hit in hits:
                        rec_a_hits_in_b[hit.data] = rec_b_map[hit.data]
                else:
                    if feature.strand > 0:
                        hits = tree_b_pos[
                            int(feature.location.start) : int(feature.location.end)
                        ]
                    else:
                        hits = tree_b_neg[
                            int(feature.location.start) : int(feature.location.end)
                        ]
                    for hit in hits:
                        rec_a_hits_in_b[hit.data] = rec_b_map[hit.data]

            for feature in rec_b_i.features:
                if feature.type == "remark" or feature.type == "annotation":
                    continue

                if not stranding:
                    hits = tree_a[
                        int(feature.location.start) : int(feature.location.end)
                    ]
                    for hit in hits:
                        rec_b_hits_in_a[hit.data] = rec_a_map[hit.data]
                else:
                    if feature.strand > 0:
                        hits = tree_a_pos[
                            int(feature.location.start) : int(feature.location.end)
                        ]
                    else:
                        hits = tree_a_neg[
                            int(feature.location.start) : int(feature.location.end)
                        ]
                    for hit in hits:
                        rec_b_hits_in_a[hit.data] = rec_a_map[hit.data]

            # Sort features by start position
            rec_a_out.append(
                SeqRecord(
                    rec_a[iterate].seq,
                    rec_a[iterate].id,
                    rec_a[iterate].name,
                    rec_a[iterate].description,
                    rec_a[iterate].dbxrefs,
                    sorted(
                        rec_a_hits_in_b.values(), key=lambda feat: feat.location.start
                    ),
                    rec_a[iterate].annotations,
                )
            )
            rec_b_out.append(
                SeqRecord(
                    rec_b[iterate].seq,
                    rec_b[iterate].id,
                    rec_b[iterate].name,
                    rec_b[iterate].description,
                    rec_b[iterate].dbxrefs,
                    sorted(
                        rec_b_hits_in_a.values(), key=lambda feat: feat.location.start
                    ),
                    rec_b[iterate].annotations,
                )
            )
            iterate += 1