graph

navis.graph.classify_nodes

Classify neuron's nodes into end nodes, branches, slabs or root.

Adds a 'type' column to x.nodes table.

PARAMETER	DESCRIPTION
x	Neuron(s) whose nodes to classify. TYPE: TreeNeuron | NeuronList
categorical	If True (default), will use categorical data type which takes up much less memory at a small run-time overhead. TYPE: bool DEFAULT: True
inplace	If `False`, nodes will be classified on a copy which is then returned leaving the original neuron unchanged. TYPE: bool DEFAULT: True

RETURNS	DESCRIPTION
TreeNeuron / List

Examples:

>>> import navis
>>> nl = navis.example_neurons(2)
>>> _ = navis.graph.classify_nodes(nl, inplace=True)

Source code in navis/graph/graph_utils.py

@utils.map_neuronlist(desc="Classifying", allow_parallel=True)
@utils.lock_neuron
def classify_nodes(x: "core.NeuronObject", categorical=True, inplace: bool = True):
    """Classify neuron's nodes into end nodes, branches, slabs or root.

    Adds a `'type'` column to `x.nodes` table.

    Parameters
    ----------
    x :         TreeNeuron | NeuronList
                Neuron(s) whose nodes to classify.
    categorical : bool
                If True (default), will use categorical data type which takes
                up much less memory at a small run-time overhead.
    inplace :   bool, optional
                If `False`, nodes will be classified on a copy which is then
                returned leaving the original neuron unchanged.

    Returns
    -------
    TreeNeuron/List

    Examples
    --------
    >>> import navis
    >>> nl = navis.example_neurons(2)
    >>> _ = navis.graph.classify_nodes(nl, inplace=True)

    """
    if not inplace:
        x = x.copy()

    if not isinstance(x, core.TreeNeuron):
        raise TypeError(f'Expected TreeNeuron(s), got "{type(x)}"')

    if x.nodes.empty:
        x.nodes["type"] = None
        return x

    node_ids = x.nodes.node_id.values
    parent_ids = x.nodes.parent_id.values

    # Note: we work with the integer *codes* of `NODE_TYPES` throughout and only
    # ever turn them into labels at the very end. Going via a string array (as
    # this used to) makes `pd.Categorical` factorize N strings, which costs more
    # than the classification itself.
    if utils.fastcore:
        # Fastcore uses its own order (0=root, 1=leaf, 2=branch, 3=slab)
        cl = _FASTCORE_NODE_TYPES[utils.fastcore.classify_nodes(node_ids, parent_ids)]
    else:
        # Note: I have tried to optimized the s**t out of this, i.e. every
        # single line of code here has been tested for speed. Do not
        # change anything unless you know what you're doing!

        # Turns out that numpy.isin() recently started to complain if the
        # node_ids are uint64 and the parent_ids are int64 (but strangely
        # not with 32bit integers). If that's the case we have to convert
        # the node_ids to int64.
        if node_ids.dtype == np.uint64:
            node_ids = node_ids.astype(np.int64)

        cl = np.full(len(x.nodes), NODE_TYPES.index("slab"), dtype=np.int8)
        cl[~np.isin(node_ids, parent_ids)] = NODE_TYPES.index("end")
        bp = x.nodes.parent_id.value_counts()
        bp = bp.index.values[bp.values > 1]
        cl[np.isin(node_ids, bp)] = NODE_TYPES.index("branch")
        cl[parent_ids < 0] = NODE_TYPES.index("root")

    if categorical:
        x.nodes["type"] = pd.Categorical.from_codes(
            cl, categories=NODE_TYPES, ordered=False
        )
    else:
        x.nodes["type"] = np.asarray(NODE_TYPES, dtype="<U6")[cl]

    return x

navis.graph.connected_components_of

Weakly connected components of the sub-graph induced on keep.

Returns sets of node IDs, mirroring nx.connected_components.

Source code in navis/graph/graph_utils.py

def connected_components_of(x: "core.TreeNeuron", keep) -> List[Set[int]]:
    """Weakly connected components of the sub-graph induced on `keep`.

    Returns sets of node IDs, mirroring `nx.connected_components`.
    """
    sub = subset_igraph(x, keep)
    ids = np.asarray(sub.vs["node_id"])
    return [set(ids[c].tolist()) for c in sub.components(mode="WEAK")]

navis.graph.connected_subgraph

Return set of nodes necessary to connect all nodes in subset ss.

PARAMETER	DESCRIPTION
x	Neuron (or graph thereof) to get subgraph for. TYPE: navis.TreeNeuron | nx.DiGraph
ss	Node IDs of node to subset to. TYPE: list | array-like

RETURNS	DESCRIPTION
np.ndarray	Node IDs of connected subgraph.
root ID	ID of the node most proximal to the old root in the connected subgraph.

Examples:

>>> import navis
>>> n = navis.example_neurons(1)
>>> ends = n.nodes[n.nodes.type.isin(['end', 'root'])].node_id.values
>>> sg, root = navis.graph.graph_utils.connected_subgraph(n, ends)
>>> # Since we asked for a subgraph connecting all terminals + root,
>>> # we expect to see all nodes in the subgraph
>>> sg.shape[0] == n.nodes.shape[0]
True

Source code in navis/graph/graph_utils.py

def connected_subgraph(
    x: Union["core.TreeNeuron", nx.DiGraph], ss: Sequence[Union[str, int]]
) -> Tuple[np.ndarray, Union[int, str]]:
    """Return set of nodes necessary to connect all nodes in subset `ss`.

    Parameters
    ----------
    x :         navis.TreeNeuron | nx.DiGraph
                Neuron (or graph thereof) to get subgraph for.
    ss :        list | array-like
                Node IDs of node to subset to.

    Returns
    -------
    np.ndarray
                Node IDs of connected subgraph.
    root ID
                ID of the node most proximal to the old root in the
                connected subgraph.

    Examples
    --------
    >>> import navis
    >>> n = navis.example_neurons(1)
    >>> ends = n.nodes[n.nodes.type.isin(['end', 'root'])].node_id.values
    >>> sg, root = navis.graph.graph_utils.connected_subgraph(n, ends)
    >>> # Since we asked for a subgraph connecting all terminals + root,
    >>> # we expect to see all nodes in the subgraph
    >>> sg.shape[0] == n.nodes.shape[0]
    True

    """
    # `src` is the TreeNeuron we can pull node/parent arrays from (if any). For a
    # bare graph we fall back to reading its edges.
    src = None
    g = None
    if isinstance(x, core.NeuronList):
        if len(x) == 1:
            src = x[0]
    elif isinstance(x, core.TreeNeuron):
        src = x
    elif isinstance(x, (nx.DiGraph, igraph.Graph)):
        g = x
    else:
        raise TypeError(f'Input must be a single TreeNeuron or graph, got "{type(x)}".')

    # Build a `node -> parent` map and the set of nodes in the graph.
    # `parent.get(n)` returns None for roots (which have no parent) - this is our
    # natural walk terminator (mirrors `next(g.successors(n), None)`).
    # For a TreeNeuron we can build this straight from the node table (faster and
    # avoids touching a graph library at all); for a bare graph (e.g. the induced
    # sub-graph passed by `split_axon_dendrite`) we do a single pass over its edges.
    if src is not None:
        nid = src.nodes.node_id.values
        pid = src.nodes.parent_id.values
        parent = {n: p for n, p in zip(nid, pid) if p >= 0}
        nodes = set(nid.tolist())
    elif isinstance(g, igraph.Graph):
        ids = np.asarray(g.vs["node_id"])
        edges = np.asarray(g.get_edgelist(), dtype=np.int64).reshape(-1, 2)
        # edge (u, v) => v is parent of u
        parent = dict(zip(ids[edges[:, 0]].tolist(), ids[edges[:, 1]].tolist()))
        nodes = set(ids.tolist())
    else:
        parent = {u: v for u, v in g.edges()}  # edge (u, v) => v is parent of u
        nodes = set(g.nodes())

    ss = set(ss)
    missing = ss - nodes
    if missing:
        missing = np.array(list(missing)).astype(str)  # do NOT remove list() here!
        raise ValueError(f"Nodes not found: {','.join(missing)}")

    # Find nodes that are leafs WITHIN the subset: an ss node is an ss-leaf iff none
    # of its children are in ss, i.e. it is not the parent of any other ss node.
    ss_parents = {parent[n] for n in ss if parent.get(n) in ss}
    leafs = ss - ss_parents

    # Memoised depth (distance to root; root = 0). Each node is resolved exactly
    # once thanks to the `n in depth` early stop -> O(N). Replaces the old
    # `longest_path.index(...)` ordering key (which was O(depth) per lookup).
    depth = {}

    def fill_depth(n):
        stack = []
        while n is not None and n not in depth:
            stack.append(n)
            n = parent.get(n)
        d = depth[n] + 1 if n is not None else 0
        for m in reversed(stack):
            depth[m] = d
            d += 1

    # Walk every ss-leaf towards its root, stopping as soon as we hit an already
    # visited node. We accumulate, per node, how many leaf-walks pass through it
    # (`pass_count`) and which component (terminal root) it belongs to (`comp_of`).
    # Components are derived implicitly: leaves ending at the same root share one.
    pass_count = {}
    comp_of = {}
    comp_leaves = defaultdict(list)
    comp_touched = defaultdict(list)
    for leaf in leafs:
        fill_depth(leaf)
        # First pass: walk to root, counting passes and finding the component root.
        n = leaf
        root = leaf
        while n is not None:
            pass_count[n] = pass_count.get(n, 0) + 1
            root = n
            n = parent.get(n)
        comp_leaves[root].append(leaf)
        # Second pass: tag every (not yet tagged) node on this path with its
        # component root and record it as touched. Early-stops where a previous
        # leaf-walk already tagged the shared upper segment.
        n = leaf
        while n is not None and n not in comp_of:
            comp_of[n] = root
            comp_touched[root].append(n)
            n = parent.get(n)

    # Group ss nodes by component once (every ss node lies on some leaf-walk and is
    # therefore tagged in `comp_of`). Avoids re-scanning all of ss per component.
    ss_by_comp = defaultdict(list)
    for n in ss:
        ss_by_comp[comp_of[n]].append(n)

    include = set()
    new_roots = []
    for root, cleaves in comp_leaves.items():
        need = len(cleaves)
        # Nodes common to ALL leaf-walks form a contiguous root->LCA chain; the LCA
        # (branch point / new root) is the deepest of them.
        common = [n for n in comp_touched[root] if pass_count[n] == need]
        lca = max(common, key=lambda n: depth[n])

        # Include, for each leaf, every node up to and including the LCA.
        for leaf in cleaves:
            n = leaf
            while n is not None and n not in include:
                include.add(n)
                if n == lca:
                    break
                n = parent.get(n)

        # Edge case: ss may contain nodes that are strict ancestors of the LCA
        # (they are never ss-leaves, so they're not in `include` yet). The new root
        # must be the most proximal of those (closest to the old root, i.e. smallest
        # depth); we then fill the *full* chain from the LCA up to it so the result
        # stays connected (the old code added only the ss nodes, leaving a
        # disconnected gap between the LCA and the new root).
        this_ss = ss_by_comp[root]
        proximal = [n for n in this_ss if n not in include]
        if proximal:
            new_root = min(proximal, key=lambda n: depth[n])
            new_roots.append(new_root)
            # All proximal ss nodes are ancestors of the LCA, so walking from the
            # LCA towards the root reaches `new_root` and passes every one of them.
            n = lca
            while True:
                include.add(n)
                if n == new_root:
                    break
                n = parent.get(n)
        else:
            new_roots.append(lca)

    return np.array(list(include)), new_roots

navis.graph.generate_list_of_childs

Return list of childs.

PARAMETER	DESCRIPTION
x	If List, must contain a SINGLE neuron. TYPE: TreeNeuron | NeuronList

RETURNS	DESCRIPTION
dict	{parent_id: [child_id, child_id, ...]}

Source code in navis/graph/graph_utils.py

def generate_list_of_childs(x: "core.NeuronObject") -> Dict[int, List[int]]:
    """Return list of childs.

    Parameters
    ----------
    x :     TreeNeuron | NeuronList
            If List, must contain a SINGLE neuron.

    Returns
    -------
    dict
        `{parent_id: [child_id, child_id, ...]}`

    """
    assert isinstance(x, core.TreeNeuron)

    # The node table already *is* a child->parent map, so we can invert it directly
    # instead of building a graph and asking it for `in_edges` once per node.
    nid = x.nodes.node_id.values
    pid = x.nodes.parent_id.values

    childs: Dict[int, List[int]] = {n: [] for n in nid.tolist()}
    has_parent = pid >= 0
    for c, p in zip(nid[has_parent].tolist(), pid[has_parent].tolist()):
        childs[p].append(c)

    return childs

navis.graph.node_label_sorting

Return nodes ordered by node label sorting according to Cuntz et al., PLoS Computational Biology (2010).

PARAMETER	DESCRIPTION
x	TYPE: TreeNeuron
weighted	If True will use actual distances instead of just node count. Depending on how evenly spaced your points are, this might not make much difference. TYPE: bool DEFAULT: False

RETURNS	DESCRIPTION
list	[root, node_id, node_id, ...]

Source code in navis/graph/graph_utils.py

def node_label_sorting(
    x: "core.TreeNeuron", weighted: bool = False
) -> List[Union[str, int]]:
    """Return nodes ordered by node label sorting according to Cuntz
    et al., PLoS Computational Biology (2010).

    Parameters
    ----------
    x :         TreeNeuron
    weighted :  bool
                If True will use actual distances instead of just node count.
                Depending on how evenly spaced your points are, this might not
                make much difference.

    Returns
    -------
    list
        `[root, node_id, node_id, ...]`

    """
    if isinstance(x, core.NeuronList) and len(x) == 1:
        x = x[0]

    if not isinstance(x, core.TreeNeuron):
        raise TypeError(f'Expected a singleTreeNeuron, got "{type(x)}"')

    if len(x.root) > 1:
        raise ValueError("Unable to process multi-root neurons!")

    # Get relevant terminal nodes
    term = x.nodes[x.nodes.type == "end"].node_id.values

    weight = "weight" if weighted else None

    # The walk below sorts each node `n` by "distance to the farthest terminal below
    # `n`" plus "distance from `n` up to the node we are walking from". This used to
    # come out of a directed breaks-by-breaks geodesic matrix, but neither term needs
    # one: the first is `n`'s subtree height and the second - since we only ever walk
    # to an ancestor of `n` - is a difference of root distances. Both are O(N).
    # The matrix was the single largest allocation in navis (19512 x 19512, i.e.
    # 4.5GB, for a 71k node skeleton).
    # Plain dicts: the sort keys below do one scalar lookup per node and pandas
    # charges ~2.5us for each of those.
    height = morpho.manipulation._subtree_height(x, weight=weight).to_dict()
    depth = dist_to_root(x, weight=weight)

    def sort_key(parent):
        """Sort a node's children by the longest path running through them."""
        return lambda n: height[n] + (depth[n] - depth[parent])

    # Get starting points (i.e. branches off the root) and sort by longest
    # path to a terminal (note we're operating on the simplified version
    # of the skeleton)
    G = graph.simplify_graph(x.graph)
    curr_points = sorted(
        list(G.predecessors(x.root[0])), key=sort_key(x.root[0]), reverse=True
    )

    # Walk from root towards terminals, prioritising longer branches
    nodes_walked = []
    while curr_points:
        nodes_walked.append(curr_points.pop(0))
        # If the current point is a terminal point, stop here
        if nodes_walked[-1] in term:
            pass
        else:
            new_points = sorted(
                list(G.predecessors(nodes_walked[-1])),
                key=sort_key(nodes_walked[-1]),
                reverse=True,
            )
            curr_points = new_points + curr_points

    # Translate into segments
    node_list = [x.root[0:]]
    # Note that we're inverting here so that the segments are ordered
    # proximal -> distal (i.e. root to tips)
    seg_dict = {s[0]: s[::-1] for s in _break_segments(x)}

    for n in nodes_walked:
        # Note that we're skipping the first (proximal) node to avoid double
        # counting nodes
        node_list.append(seg_dict[n][1:])

    return np.concatenate(node_list, dtype=int)

navis.graph.segment_lengths

Get the length of each of many linear segments.

Same as calling navis.segment_length on each segment but builds the node lookup once instead of once per segment.

RETURNS	DESCRIPTION
np.ndarray	Length of each segment.

Source code in navis/graph/graph_utils.py

def segment_lengths(x: "core.TreeNeuron", segments: Sequence[Sequence[int]]):
    """Get the length of each of many linear segments.

    Same as calling [`navis.segment_length`][] on each segment but builds the node
    lookup once instead of once per segment.

    Returns
    -------
    np.ndarray
                Length of each segment.
    """
    if not len(segments):
        return np.zeros(0)

    # An edge's weight is just the distance between its two nodes, so we can read
    # the lengths straight off the coordinates rather than going via a graph.
    # Note the cast to float64: node coordinates are often float32, and summing
    # those would drift away from the weights networkx used to hand us.
    coords = x.nodes[["x", "y", "z"]].values.astype(float)

    # Resolve every segment's node IDs in one lookup - `get_indexer` has enough
    # per-call overhead that doing it once per segment costs more than the walk it
    # replaces.
    lengths = np.array([len(s) for s in segments])
    flat = np.concatenate([np.asarray(s) for s in segments])
    coords = coords[pd.Index(x.nodes.node_id).get_indexer(flat)]

    # Distance from each node to the one before it...
    step = np.zeros(len(flat))
    step[1:] = np.linalg.norm(np.diff(coords, axis=0), axis=1)

    # ...except the first node of each segment, which has no predecessor *in that
    # segment* - this also discards the bogus step across each segment boundary.
    starts = np.concatenate([[0], np.cumsum(lengths)[:-1]])
    step[starts] = 0

    return np.add.reduceat(step, starts)

navis.graph.simplify_graph

Simplify skeleton graph (networkX or igraph).

This function will simplify the graph by keeping only roots, leafs and branch points. Preserves branch lengths (i.e. weights)!

PARAMETER	DESCRIPTION
G	The skeleton graph to simplify. TYPE: networkx.DiGraph | igraph.Graph
inplace	If True, will modify the graph in place. TYPE: bool DEFAULT: False

RETURNS	DESCRIPTION
G	Simplified graph. TYPE: networkx.DiGraph | networkx.DiGraph

Examples:

>>> import navis
>>> n = navis.example_neurons(1, kind='skeleton')
>>> # Simplify skeleton's NetworkX graph representation
>>> G_simp_nx = navis.graph.simplify_graph(n.graph)
>>> # Check that we have the expected number of nodes
>>> assert len(G_simp_nx.nodes) == (n.n_branches + n.n_root + n.n_leafs)
>>> # Simplify skeleton's iGraph graph representation
>>> G_simp_ig = navis.graph.simplify_graph(n.igraph)
>>> # Check that we have the expected number of nodes
>>> assert len(G_simp_ig.vs) == (n.n_branches + n.n_root + n.n_leafs)

Source code in navis/graph/converters.py

def simplify_graph(G, inplace=False):
    """Simplify skeleton graph (networkX or igraph).

    This function will simplify the graph by keeping only roots, leafs and
    branch points. Preserves branch lengths (i.e. weights)!

    Parameters
    ----------
    G :         networkx.DiGraph | igraph.Graph
                The skeleton graph to simplify.
    inplace :   bool
                If True, will modify the graph in place.

    Returns
    -------
    G :         networkx.DiGraph | networkx.DiGraph
                Simplified graph.

    Examples
    --------
    >>> import navis
    >>> n = navis.example_neurons(1, kind='skeleton')
    >>> # Simplify skeleton's NetworkX graph representation
    >>> G_simp_nx = navis.graph.simplify_graph(n.graph)
    >>> # Check that we have the expected number of nodes
    >>> assert len(G_simp_nx.nodes) == (n.n_branches + n.n_root + n.n_leafs)
    >>> # Simplify skeleton's iGraph graph representation
    >>> G_simp_ig = navis.graph.simplify_graph(n.igraph)
    >>> # Check that we have the expected number of nodes
    >>> assert len(G_simp_ig.vs) == (n.n_branches + n.n_root + n.n_leafs)

    """
    if not inplace:
        G = G.copy()

    if isinstance(G, nx.Graph):
        # Find all leaf and branch points
        leafs = {n for n in G.nodes if G.in_degree(n) == 0 and G.out_degree(n) != 0}
        branches = {n for n in G.nodes if G.in_degree(n) > 1 and G.out_degree(n) != 0}
        roots = {n for n in G.nodes if G.out_degree(n) == 0}

        stop_nodes = roots | leafs | branches

        # Walk from each leaf/branch point to the next leaf, branch or root
        to_remove = []
        for start_node in leafs | branches:
            dist = 0
            node = start_node
            while True:
                parent = next(G.successors(node))
                dist += G.edges[node, parent]["weight"]

                if parent in stop_nodes:
                    G.add_weighted_edges_from([(start_node, parent, dist)])
                    break

                to_remove.append(parent)
                node = parent

        G.remove_nodes_from(to_remove)
    else:
        # Find all leaf and branch points
        leafs = G.vs.select(_indegree=0, _outdegree_ne=0)
        branches = G.vs.select(_indegree_gt=1, _outdegree_ne=0)
        roots = G.vs.select(_outdegree=0)

        stop_nodes = np.concatenate((roots.indices, leafs.indices, branches.indices))

        # Walk from each leaf/branch point to the next leaf, branch or root
        to_remove = []
        for start_node in np.concatenate((leafs.indices, branches.indices)):
            dist = 0
            node = start_node
            while True:
                parent = G.successors(node)[0]
                dist += G.es[G.get_eid(node, parent)]["weight"]

                if parent in stop_nodes:
                    G.add_edge(start_node, parent, weight=dist)
                    break

                to_remove.append(parent)
                node = parent

        G.delete_vertices(to_remove)

    return G

navis.graph.skeleton_adjacency_matrix

Generate adjacency matrix for a skeleton.

PARAMETER	DESCRIPTION
x	Neuron for which to generate adjacency matrix. TYPE: TreeNeuron
sort	If True, will sort the adjacency matrix by topology. TYPE: bool DEFAULT: True

RETURNS	DESCRIPTION
pd.DataFrame	Adjacency matrix where rows are nodes and columns are their parents.

See Also

navis.geodesic_matrix For distances between all points. navis.distal_to Check if a node A is distal to node B. navis.dist_between Get point-to-point geodesic ("along-the-arbor") distances.

Source code in navis/graph/graph_utils.py

def skeleton_adjacency_matrix(
    x: "core.NeuronObject", sort: bool = True
) -> pd.DataFrame:
    """Generate adjacency matrix for a skeleton.

    Parameters
    ----------
    x :         TreeNeuron
                Neuron for which to generate adjacency matrix.
    sort :      bool, optional
                If True, will sort the adjacency matrix by topology.

    Returns
    -------
    pd.DataFrame
                Adjacency matrix where rows are nodes and columns are
                their parents.

    See Also
    --------
    [`navis.geodesic_matrix`][]
        For distances between all points.
    [`navis.distal_to`][]
        Check if a node A is distal to node B.
    [`navis.dist_between`][]
        Get point-to-point geodesic ("along-the-arbor") distances.

    """
    if isinstance(x, core.NeuronList):
        if len(x) == 1:
            x = x[0]
        else:
            raise ValueError("Cannot process more than a single neuron.")
    elif not isinstance(x, (core.TreeNeuron,)):
        raise ValueError(f'Unable to process data of type "{type(x)}"')

    # Generate the empty adjacency matrix
    adj = pd.DataFrame(
        np.zeros((len(x.nodes), len(x.nodes)), dtype=bool),
        index=x.nodes.node_id.values,
        columns=x.nodes.node_id.values,
    )

    # Fill in the parent-child relationships
    not_root = x.nodes.parent_id.values >= 0
    node_ix = np.arange(len(x.nodes))[not_root]
    parent_ids = x.nodes.parent_id.values[not_root]
    parent_ix = np.searchsorted(x.nodes.node_id.values, parent_ids)
    adj.values[node_ix, parent_ix] = True

    if sort:
        sort = node_label_sorting(x)
        adj = adj.loc[sort, sort]

    return adj

navis.graph.subset_igraph

Induce the sub-graph of a neuron's igraph on a set of node IDs.

The igraph equivalent of x.graph.subgraph(keep), but without paying to build the networkx graph. node_id is carried over onto the new vertices.

Source code in navis/graph/graph_utils.py

def subset_igraph(x: "core.TreeNeuron", keep) -> "igraph.Graph":
    """Induce the sub-graph of a neuron's igraph on a set of node IDs.

    The igraph equivalent of `x.graph.subgraph(keep)`, but without paying to build
    the networkx graph. `node_id` is carried over onto the new vertices.
    """
    G: igraph.Graph = x.igraph
    ids = np.asarray(G.vs["node_id"])
    keep = np.fromiter(keep, dtype=ids.dtype, count=len(keep))
    return G.subgraph(np.where(np.isin(ids, keep))[0])