From 4c774a3ab24c79f2d1d87925cd76b7c4ff031f98 Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Wed, 12 Aug 2026 16:25:31 -0400
Subject: [PATCH 01/11] MultiOrderModel now holds HigherOrderGraph in its
 layers

---
 src/pathpyG/__init__.py               |  2 +
 src/pathpyG/core/multi_order_model.py | 69 +++++++++++++++++++--------
 2 files changed, 50 insertions(+), 21 deletions(-)

diff --git a/src/pathpyG/__init__.py b/src/pathpyG/__init__.py
index 4c365455..c4c23c9a 100644
--- a/src/pathpyG/__init__.py
+++ b/src/pathpyG/__init__.py
@@ -8,6 +8,7 @@
 __version__ = get_version("pathpyG")
 
 from pathpyG.core.graph import Graph
+from pathpyG.core.higher_order_graph import HigherOrderGraph
 from pathpyG.core.index_map import IndexMap
 from pathpyG.core.multi_order_model import MultiOrderModel
 from pathpyG.core.path_data import PathData
@@ -21,6 +22,7 @@
 
 __all__ = [
     "Graph",
+    "HigherOrderGraph",
     "TemporalGraph",
     "EventGraph",
     "PathData",
diff --git a/src/pathpyG/core/multi_order_model.py b/src/pathpyG/core/multi_order_model.py
index 542156ac..88d068f7 100644
--- a/src/pathpyG/core/multi_order_model.py
+++ b/src/pathpyG/core/multi_order_model.py
@@ -17,11 +17,10 @@
     lift_order_edge_index_weighted,
 )
 from pathpyG.core.event_graph import EventGraph
-from pathpyG.core.graph import Graph
+from pathpyG.core.higher_order_graph import HigherOrderGraph
 from pathpyG.core.index_map import IndexMap
 from pathpyG.core.path_data import PathData
 from pathpyG.core.temporal_graph import TemporalGraph
-from pathpyG.utils.dbgnn import generate_bipartite_edge_index
 
 logger = logging.getLogger("root")
 
@@ -31,13 +30,17 @@ class MultiOrderModel:
 
     This class stores multiple higher-order De Bruijn graphs as layers in a dictionary.
     Each layer corresponds to a De Bruijn graph of order k, where k is the key in the dictionary.
-    Each graph layer is represented as a [pathpyG.Graph][] object.
+    Each graph layer is represented as a
+    [HigherOrderGraph][pathpyG.core.higher_order_graph.HigherOrderGraph] object, layer 1
+    included. Each layer therefore knows its own order and the first-order nodes it was
+    built from, so results can be projected back onto entities without the caller keeping
+    track of it.
     This class provides methods to search for the optimal order of the model based on likelihood ratio tests,
     as well as methods to compute the log-likelihood of observed paths given the model.
 
     Attributes:
-        layers (dict[int, Graph]): A dictionary mapping the order k to the corresponding
-            higher-order De Bruijn graph of order k.
+        layers (dict[int, HigherOrderGraph]): A dictionary mapping the order k to the
+            corresponding higher-order De Bruijn graph of order k.
 
     Examples:
         Example where the optimal order is 1:
@@ -56,11 +59,15 @@ class MultiOrderModel:
         >>> m = MultiOrderModel.from_path_data(paths, max_order=2)
         >>> print(m.estimate_order(paths, max_order=2))
         2
+
+        Each layer knows the first-order path that each of its nodes represents:
+        >>> print(m.layers[2].order, m.layers[2].nodes)
+        2 [('a', 'c'), ('b', 'c'), ('c', 'd'), ('c', 'e')]
     """
 
     def __init__(self) -> None:
         """Initialize an empty MultiOrderModel."""
-        self.layers: dict[int, Graph] = {}
+        self.layers: dict[int, HigherOrderGraph] = {}
 
     def __str__(self) -> str:
         """Return a string representation of the higher-order graph."""
@@ -88,7 +95,8 @@ def iterate_lift_order(
         edge_weight: torch.Tensor | None = None,
         aggr: str = "src",
         save: bool = True,
-    ) -> tuple[torch.Tensor, torch.Tensor, torch.Tensor | None, Graph | None]:
+        n_first_order: Optional[int] = None,
+    ) -> tuple[torch.Tensor, torch.Tensor, torch.Tensor | None, HigherOrderGraph | None]:
         """Lift order by one and save the result in the layers dictionary of the object.
 
         This is a helper function that should not be called directly.
@@ -103,6 +111,8 @@ def iterate_lift_order(
             k: The order of the graph that should be computed.
             aggr: The aggregation method to use. One of "src", "dst", "max", "mul".
             save: Whether to compute the aggregated graph and later save it in the layers dictionary.
+            n_first_order: The number of first-order nodes the node sequences refer to.
+                Defaults to the number of IDs in `mapping`.
         """
         # Lift order
         if edge_weight is None:
@@ -115,8 +125,13 @@ def iterate_lift_order(
 
         # Aggregate
         if save:
-            gk = aggregate_edge_index(ho_index, node_sequence, edge_weight)
-            gk.mapping = IndexMap([tuple(mapping.to_ids(v.cpu())) for v in gk.data.node_sequence])
+            gk = HigherOrderGraph.from_aggregated(
+                ho_index,
+                node_sequence,
+                first_order_mapping=mapping,
+                edge_weight=edge_weight,
+                n_first_order=n_first_order,
+            )
         else:
             gk = None
         return ho_index, node_sequence, edge_weight, gk
@@ -156,10 +171,13 @@ def from_temporal_graph(
         else:
             edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
         if cached or max_order == 1:
-            m.layers[1] = aggregate_edge_index(
-                edge_index=edge_index, node_sequence=node_sequence, edge_weight=edge_weight
+            m.layers[1] = HigherOrderGraph.from_aggregated(
+                edge_index=edge_index,
+                node_sequence=node_sequence,
+                edge_weight=edge_weight,
+                first_order_mapping=g.mapping,
+                n_first_order=g.n,
             )
-            m.layers[1].mapping = g.mapping
 
         if max_order > 1:
             node_sequence = torch.cat([node_sequence[edge_index[0]], node_sequence[edge_index[1]][:, -1:]], dim=1)
@@ -171,11 +189,12 @@ def from_temporal_graph(
 
             # Aggregate
             if cached or max_order == 2:
-                m.layers[2] = aggregate_edge_index(
-                    edge_index=edge_index, node_sequence=node_sequence, edge_weight=edge_weight
-                )
-                m.layers[2].mapping = IndexMap(
-                    [tuple(g.mapping.to_ids(v.cpu())) for v in m.layers[2].data.node_sequence]
+                m.layers[2] = HigherOrderGraph.from_aggregated(
+                    edge_index=edge_index,
+                    node_sequence=node_sequence,
+                    edge_weight=edge_weight,
+                    first_order_mapping=g.mapping,
+                    n_first_order=g.n,
                 )
 
             for k in range(3, max_order + 1):
@@ -186,6 +205,7 @@ def from_temporal_graph(
                     edge_weight=edge_weight,
                     aggr="src",
                     save=cached or k == max_order,
+                    n_first_order=g.n,
                 )
                 if cached or k == max_order:
                     m.layers[k] = gk  # type: ignore[assignment]
@@ -251,17 +271,24 @@ def from_path_data(
         elif mode == "propagation":
             aggr = "src"
 
-        m.layers[1] = aggregate_edge_index(edge_index=edge_index, node_sequence=node_sequence, edge_weight=edge_weight)
-        m.layers[1].mapping = path_data.mapping
+        g1 = aggregate_edge_index(edge_index=edge_index, node_sequence=node_sequence, edge_weight=edge_weight)
+        g1.mapping = path_data.mapping
+        # Nodes that are not traversed by any path are not part of the aggregated graph,
+        # so the first-order node set can be larger than the order-1 layer.
+        n_first_order = max(path_data.mapping.num_ids(), g1.n)
+        m.layers[1] = HigherOrderGraph.from_aggregated_graph(
+            g1, first_order_mapping=path_data.mapping, n_first_order=n_first_order
+        )
 
         for k in range(2, max_order + 1):
             edge_index, node_sequence, edge_weight, gk = MultiOrderModel.iterate_lift_order(
                 edge_index=edge_index,
                 node_sequence=node_sequence,
-                mapping=m.layers[1].mapping,
+                mapping=path_data.mapping,
                 edge_weight=edge_weight,
                 aggr=aggr,
                 save=cached or k == max_order,
+                n_first_order=n_first_order,
             )
             if cached or k == max_order:
                 m.layers[k] = gk  # type: ignore[assignment]
@@ -563,7 +590,7 @@ def to_dbgnn_data(self, max_order: int = 2, mapping: str = "last") -> Data:
         edge_index_max_order = g_max_order.data.edge_index
         edge_weight = g.data.edge_weight
         edge_weight_max_order = g_max_order.data.edge_weight
-        bipartite_edge_index = generate_bipartite_edge_index(g, g_max_order, mapping=mapping, device=edge_index.device)
+        bipartite_edge_index = g_max_order.bipartite_edge_index(g, mapping=mapping, device=edge_index.device)
 
         if g.data.y is not None:
             y = g.data.y

From 59bfa8b695f9adc587ce180810ff96077297fc1d Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Wed, 12 Aug 2026 16:59:02 -0400
Subject: [PATCH 02/11] added missing HigherOrderGraph class

---
 src/pathpyG/core/higher_order_graph.py | 521 +++++++++++++++++++++++++
 1 file changed, 521 insertions(+)
 create mode 100644 src/pathpyG/core/higher_order_graph.py

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
new file mode 100644
index 00000000..56b7bebc
--- /dev/null
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -0,0 +1,521 @@
+"""Higher-order De Bruijn graph representation and related operations."""
+
+from __future__ import annotations
+
+import logging
+from typing import Optional, Union
+
+import torch
+from torch_geometric import EdgeIndex
+from torch_geometric.data import Data
+from torch_geometric.utils import coalesce
+
+from pathpyG.algorithms.lift_order import aggregate_edge_index, lift_order_edge_index_weighted
+from pathpyG.core.event_graph import EventGraph
+from pathpyG.core.graph import Graph
+from pathpyG.core.index_map import IndexMap
+from pathpyG.core.path_data import PathData
+from pathpyG.core.temporal_graph import TemporalGraph
+
+logger = logging.getLogger("root")
+
+
+class HigherOrderGraph(Graph):
+    """A De Bruijn graph of order `k`, whose nodes are paths of `k` first-order nodes.
+
+    Where a [`Graph`][pathpyG.Graph] has one node per entity and an
+    [`EventGraph`][pathpyG.core.event_graph.EventGraph] has one node per observed
+    interaction, a `HigherOrderGraph` has one node per *distinct* path of length `k`
+    in the underlying first-order graph. Repeated observations of the same path are
+    aggregated into an `edge_weight`, so timestamps are no longer represented: this
+    is a model of how paths flow rather than a record of what happened.
+
+    Order 1 is the degenerate case and is simply the weighted first-order graph, with
+    plain node IDs rather than tuples.
+
+    Info:
+        In addition to the attributes of [`Graph`][pathpyG.Graph], the `data` object holds:
+
+        - `node_sequence`: [Tensor][torch.Tensor] of shape `(num_nodes, order)`, the
+            first-order node indices making up the path each higher-order node represents.
+        - `edge_weight`: [Tensor][torch.Tensor] with the aggregated weight of each transition.
+        - `inverse_idx`: [Tensor][torch.Tensor] mapping each row of the *pre-aggregation*
+            node sequence to the index of the higher-order node it was merged into.
+
+    Attributes:
+        data (Data): PyG Data object containing edges and attributes.
+        mapping (IndexMap): Mapping from higher-order node IDs (tuples, for order > 1) to indices.
+        first_order_mapping (IndexMap): Mapping of the underlying first-order node IDs to indices.
+        n_first_order (int): Number of first-order nodes the higher-order nodes are built from.
+
+    Examples:
+        >>> import pathpyG as pp
+        >>> from pathpyG.core.higher_order_graph import HigherOrderGraph
+        >>> g = pp.Graph.from_edge_list([("a", "c"), ("c", "d")])
+        >>> h = HigherOrderGraph.from_graph(g)
+        >>> print(h.order, h.nodes)
+        1 ['a', 'c', 'd']
+    """
+
+    def __init__(
+        self,
+        data: Data,
+        order: Optional[int] = None,
+        first_order_mapping: Optional[IndexMap] = None,
+        n_first_order: Optional[int] = None,
+        mapping: Optional[IndexMap] = None,
+    ) -> None:
+        """Create a HigherOrderGraph from a `Data` object carrying a `node_sequence`.
+
+        Args:
+            data: PyG `Data` object with an `edge_index` and a `node_sequence` of shape
+                `(num_nodes, order)`. For order 1, the `node_sequence` may be omitted and
+                is then taken to be the identity.
+            order: Expected order `k`. If given, it is validated against the width of the
+                node sequence; if omitted, the order is inferred from it.
+            first_order_mapping: Mapping of the underlying first-order node IDs. Defaults
+                to an empty mapping.
+            n_first_order: Number of first-order nodes. Defaults to the number of IDs in
+                `first_order_mapping`, or the largest index in the node sequence plus one.
+            mapping: Mapping of higher-order node IDs to indices. For order > 1 this must
+                use tuple IDs; for order 1 it must not.
+
+        Raises:
+            ValueError: If the order, the node sequence, and the mapping disagree, or if
+                the node sequence refers to first-order nodes that do not exist.
+        """
+        if "node_sequence" not in data and order not in (None, 1):
+            raise ValueError(f"A HigherOrderGraph of order {order} requires a `node_sequence` node attribute.")
+
+        if isinstance(data.edge_index, EdgeIndex):
+            # `Graph.__init__` re-sorts the edge index and reindexes every edge attribute by
+            # the returned permutation - but `EdgeIndex.sort_by` returns `None` for an index
+            # already known to be sorted, and `attr[None]` would add a dimension. Higher-order
+            # graphs are routinely built from already-aggregated (hence sorted) data, so hand
+            # the base class a plain tensor and let it derive a real permutation.
+            data.edge_index = data.edge_index.as_tensor()
+
+        super().__init__(data, mapping=mapping)
+
+        # `Graph` creates an identity node sequence if none is given, so `self.order`
+        # (inherited: the width of the node sequence) is now well-defined.
+        if order is not None and order != self.order:
+            raise ValueError(f"order={order} does not match node sequence of width {self.order}")
+
+        if first_order_mapping is not None:
+            self.first_order_mapping = first_order_mapping
+        elif self.order == 1:
+            # For order 1 the higher-order nodes *are* the first-order nodes.
+            self.first_order_mapping = self.mapping
+        else:
+            self.first_order_mapping = IndexMap()
+
+        if n_first_order is not None:
+            self._n_first_order = int(n_first_order)
+        elif self.first_order_mapping.has_ids:
+            self._n_first_order = self.first_order_mapping.num_ids()
+        elif self.data.node_sequence.numel() > 0:
+            self._n_first_order = int(self.data.node_sequence.max().item()) + 1
+        else:
+            self._n_first_order = 0
+
+        self._validate()
+
+    def _validate(self) -> None:
+        """Check that order, node sequence, mapping and first-order node set agree."""
+        if self.data.node_sequence.numel() > 0:
+            max_idx = int(self.data.node_sequence.max().item())
+            if max_idx >= self._n_first_order:
+                raise ValueError(
+                    f"node sequence refers to first-order node {max_idx}, "
+                    f"but there are only {self._n_first_order} first-order nodes"
+                )
+
+        if self.mapping.has_ids:
+            # Higher-order nodes are paths and are identified by tuples; first-order
+            # nodes are entities and are identified by plain IDs.
+            if self.mapping.has_tuple_ids != (self.order > 1):
+                raise ValueError(
+                    f"a mapping for a graph of order {self.order} must "
+                    f"{'use' if self.order > 1 else 'not use'} tuple IDs"
+                )
+            if self.mapping.num_ids() != self.n:
+                logger.warning(
+                    "mapping has %s IDs but graph has %s nodes", self.mapping.num_ids(), self.n
+                )
+
+    @staticmethod
+    def _validate_order(order: int) -> None:
+        """Reject orders for which no De Bruijn graph is defined."""
+        if order < 1:
+            logger.error("order must be at least 1, got %s", order)
+            raise ValueError(f"order must be at least 1, got {order}")
+
+    @staticmethod
+    def _build_mapping(node_sequence: torch.Tensor, first_order_mapping: IndexMap) -> IndexMap:
+        """Build the higher-order `IndexMap` naming each node by the path it represents."""
+        # TODO: Is it better to have a single HigherOrderMapping class?
+        order = node_sequence.size(1)
+        if node_sequence.size(0) == 0:
+            # An order beyond the longest observed path yields a graph without nodes,
+            # and `IndexMap` cannot be built from an empty list of IDs.
+            return IndexMap()
+        if order == 1:
+            # Order-1 node indices are first-order node indices, so the mapping carries over.
+            return first_order_mapping
+        if first_order_mapping.has_ids:
+            return IndexMap([tuple(first_order_mapping.to_ids(v.cpu())) for v in node_sequence])
+        return IndexMap([tuple(v.tolist()) for v in node_sequence])
+
+    @classmethod
+    def from_aggregated(
+        cls,
+        edge_index: torch.Tensor,
+        node_sequence: torch.Tensor,
+        first_order_mapping: Optional[IndexMap] = None,
+        edge_weight: Optional[torch.Tensor] = None,
+        n_first_order: Optional[int] = None,
+        aggr: str = "sum",
+    ) -> HigherOrderGraph:
+        """Aggregate a (possibly duplicated) higher-order edge index into a De Bruijn graph.
+
+        This is the single place where higher-order nodes get their identity: duplicate
+        node sequences are merged, edge weights are aggregated, and the higher-order
+        `IndexMap` naming each node by its path is built.
+
+        Args:
+            edge_index: Edge index whose nodes are indices into `node_sequence`.
+            node_sequence: Tensor of shape `(num_nodes, order)` with the first-order path
+                each (not yet aggregated) node represents.
+            first_order_mapping: Mapping of the underlying first-order node IDs.
+            edge_weight: Weight of each edge prior to aggregation. Defaults to ones.
+            n_first_order: Number of first-order nodes, including isolated ones.
+            aggr: Reduction used for the edge weights. One of "sum", "mean", "min", "max".
+
+        Returns:
+            HigherOrderGraph: The aggregated higher-order graph.
+        """
+        if isinstance(edge_index, torch.Tensor) and hasattr(edge_index, "as_tensor"):
+            edge_index = edge_index.as_tensor()
+
+        order = node_sequence.size(1)
+        if first_order_mapping is None:
+            first_order_mapping = IndexMap()
+        if n_first_order is None:
+            if first_order_mapping.has_ids:
+                n_first_order = first_order_mapping.num_ids()
+            else:
+                n_first_order = int(node_sequence.max().item()) + 1 if node_sequence.numel() > 0 else 0
+
+        data = aggregate_edge_index(edge_index, node_sequence, edge_weight, aggr=aggr).data
+
+        if order == 1 and n_first_order > data.num_nodes:
+            # Order-1 indices are first-order indices, so first-order nodes that are not
+            # traversed by any path are simply isolated nodes of the order-1 graph.
+            data.num_nodes = n_first_order
+            data.node_sequence = torch.arange(n_first_order, device=edge_index.device).unsqueeze(1)
+
+        return cls(
+            data,
+            order=order,
+            first_order_mapping=first_order_mapping,
+            n_first_order=n_first_order,
+            mapping=cls._build_mapping(data.node_sequence, first_order_mapping),
+        )
+
+    @classmethod
+    def from_aggregated_graph(
+        cls,
+        g: Graph,
+        first_order_mapping: Optional[IndexMap] = None,
+        n_first_order: Optional[int] = None,
+    ) -> HigherOrderGraph:
+        """Adopt an already-aggregated [`Graph`][pathpyG.Graph] as a higher-order graph.
+
+        Used to give the layers computed by a multi-order model their proper type. The
+        underlying `data` object is shared, not copied.
+
+        Args:
+            g: Aggregated graph carrying a `node_sequence` of shape `(num_nodes, order)`.
+            first_order_mapping: Mapping of the underlying first-order node IDs.
+            n_first_order: Number of first-order nodes.
+
+        Returns:
+            HigherOrderGraph: The same graph, typed as a higher-order graph.
+        """
+        if isinstance(g, HigherOrderGraph):
+            return g
+        return cls(
+            g.data,
+            first_order_mapping=first_order_mapping,
+            n_first_order=n_first_order,
+            mapping=g.mapping,
+        )
+
+    @classmethod
+    def from_graph(cls, g: Graph, weight: str = "edge_weight") -> HigherOrderGraph:
+        """Create the order-1 graph corresponding to a first-order graph.
+
+        Multi-edges are coalesced into a single weighted edge.
+
+        Args:
+            g: First-order graph.
+            weight: Name of the edge attribute to use as edge weight. If absent, each
+                edge counts once.
+
+        Returns:
+            HigherOrderGraph: A higher-order graph of order 1.
+        """
+        edge_index = g.data.edge_index.as_tensor()
+        if weight in g.data:
+            edge_weight = g.data[weight]
+        else:
+            edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
+        node_sequence = torch.arange(g.n, device=edge_index.device).unsqueeze(1)
+
+        return cls.from_aggregated(
+            edge_index,
+            node_sequence,
+            first_order_mapping=g.mapping,
+            edge_weight=edge_weight,
+            n_first_order=g.n,
+        )
+
+    @classmethod
+    def from_temporal_graph(
+        cls,
+        g: TemporalGraph,
+        order: int = 1,
+        delta: float | int = 1,
+        weight: str = "edge_weight",
+    ) -> HigherOrderGraph:
+        """Create the De Bruijn graph of order `k` for time-respecting paths in a temporal graph.
+
+        Order 1 is simply the weighted static graph and ignores `delta`; for higher orders
+        the nodes are the time-respecting paths of `k` nodes, i.e. those whose consecutive
+        interactions are at most `delta` apart. Orders above 2 are reached by repeatedly
+        lifting the *unaggregated* data, so the edge weights count observed paths rather
+        than being implied by lower-order statistics (unlike [`lift`][pathpyG.HigherOrderGraph.lift]).
+
+        Args:
+            g: The temporal graph.
+            order: The order `k` of the graph to compute.
+            delta: The maximum time difference between two consecutive interactions of a path.
+            weight: The edge attribute of `g` to use as edge weight.
+
+        Returns:
+            HigherOrderGraph: A higher-order graph of order `order`. It has no nodes if
+            there is no time-respecting path of that length.
+
+        Examples:
+            >>> import pathpyG as pp
+            >>> t = pp.TemporalGraph.from_edge_list([("a", "c", 1), ("c", "d", 2)])
+            >>> print(pp.HigherOrderGraph.from_temporal_graph(t, order=2, delta=1).nodes)
+            [('a', 'c'), ('c', 'd')]
+        """
+        cls._validate_order(order)
+        # Imported here because `MultiOrderModel` builds `HigherOrderGraph` layers itself.
+        from pathpyG.core.multi_order_model import MultiOrderModel
+
+        return MultiOrderModel.from_temporal_graph(
+            g, delta=delta, max_order=order, weight=weight, cached=False
+        ).layers[order]
+
+    @classmethod
+    def from_path_data(cls, path_data: PathData, order: int = 1, mode: str = "propagation") -> HigherOrderGraph:
+        """Create the De Bruijn graph of order `k` modelling paths in [`PathData`][pathpyG.PathData].
+
+        Args:
+            path_data: The observed paths.
+            order: The order `k` of the graph to compute.
+            mode: The process that we assume. Either "diffusion" or "propagation".
+
+        Returns:
+            HigherOrderGraph: A higher-order graph of order `order`. It has no nodes if
+            no observed path is that long.
+
+        Examples:
+            >>> import pathpyG as pp
+            >>> paths = pp.PathData(pp.IndexMap(list("acd")))
+            >>> paths.append_walk(("a", "c", "d"), weight=2)
+            >>> print(pp.HigherOrderGraph.from_path_data(paths, order=2).nodes)
+            [('a', 'c'), ('c', 'd')]
+        """
+        cls._validate_order(order)
+        from pathpyG.core.multi_order_model import MultiOrderModel
+
+        return MultiOrderModel.from_path_data(path_data, max_order=order, mode=mode, cached=False).layers[order]
+
+    @classmethod
+    def from_event_graph(cls, eg: EventGraph, order: int = 2) -> HigherOrderGraph:
+        """Aggregate an [`EventGraph`][pathpyG.core.event_graph.EventGraph] into an order-`k` graph.
+
+        For the default order 2, every event whose underlying `(u, v)` pair is the same
+        collapses into a single second-order node, and repeated continuations become an
+        edge weight. Timestamps and the time window `delta` are not represented in the
+        result. Other orders are computed from the time-respecting paths that the event
+        graph encodes, which for orders above 2 means lifting its continuations further.
+
+        Args:
+            eg: The second-order temporal event graph to aggregate.
+            order: The order `k` of the graph to compute.
+
+        Returns:
+            HigherOrderGraph: A higher-order graph of order `order`. It has no nodes if
+            there is no time-respecting path of that length.
+        """
+        if order != 2:
+            cls._validate_order(order)
+            from pathpyG.core.multi_order_model import MultiOrderModel
+
+            return MultiOrderModel.from_event_graph(eg, max_order=order, cached=False).layers[order]
+
+        edge_index = eg.data.edge_index.as_tensor()
+        # Each continuation carries the weight of the event it starts from, matching the
+        # "src" aggregation used when building order-2 layers from a temporal graph.
+        edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
+
+        return cls.from_aggregated(
+            edge_index,
+            eg.data.node_sequence,
+            first_order_mapping=eg.first_order_mapping,
+            edge_weight=edge_weight,
+            n_first_order=eg.n_first_order,
+        )
+
+    def lift(self, aggr: str = "src") -> HigherOrderGraph:
+        """Return the De Bruijn graph of order `k + 1` obtained by lifting this graph.
+
+        Nodes of the result are the edges of this graph, i.e. the paths of length `k + 1`
+        that exist in this graph's topology.
+
+        Warning:
+            This lifts an *aggregated* graph, so the resulting edge weights are those
+            implied by the order-`k` statistics rather than counts of observed paths of
+            length `k + 1`. To fit a layer to observations, use
+            [`MultiOrderModel`][pathpyG.MultiOrderModel], which lifts the unaggregated data.
+
+        Args:
+            aggr: Aggregation used for the lifted edge weights. One of "src", "dst",
+                "max", "mul" or "add".
+
+        Returns:
+            HigherOrderGraph: A higher-order graph of order `k + 1`.
+        """
+        edge_index = self.data.edge_index.as_tensor()
+        if "edge_weight" in self.data:
+            edge_weight = self.data.edge_weight
+        else:
+            edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
+
+        ho_index, ho_weight = lift_order_edge_index_weighted(
+            edge_index, edge_weight=edge_weight, num_nodes=self.n, aggr=aggr
+        )
+        node_sequence = torch.cat(
+            [self.data.node_sequence[edge_index[0]], self.data.node_sequence[edge_index[1]][:, -1:]], dim=1
+        )
+
+        return HigherOrderGraph.from_aggregated(
+            ho_index,
+            node_sequence,
+            first_order_mapping=self.first_order_mapping,
+            edge_weight=ho_weight,
+            n_first_order=self.n_first_order,
+        )
+
+    def to_first_order(self, mode: str = "last") -> Graph:
+        """Project the higher-order graph back onto the first-order nodes.
+
+        Each higher-order node is replaced by one of the first-order nodes of its path,
+        and the weights of higher-order edges mapping to the same first-order edge are
+        summed. First-order nodes not traversed by any path remain as isolated nodes.
+
+        Args:
+            mode: Which first-order node of the path represents it. Either "last" or "first".
+
+        Returns:
+            Graph: A weighted first-order graph.
+        """
+        if mode == "last":
+            projection = self.data.node_sequence[:, -1]
+        elif mode == "first":
+            projection = self.data.node_sequence[:, 0]
+        else:
+            raise ValueError(f"Unknown mode {mode}. Only 'last' and 'first' are accepted.")
+
+        edge_index = projection[self.data.edge_index.as_tensor()]
+        if "edge_weight" in self.data:
+            edge_weight = self.data.edge_weight
+        else:
+            edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
+        edge_index, edge_weight = coalesce(
+            edge_index, edge_attr=edge_weight, num_nodes=self.n_first_order, reduce="sum"
+        )
+
+        return Graph(
+            Data(edge_index=edge_index, edge_weight=edge_weight, num_nodes=self.n_first_order),
+            mapping=self.first_order_mapping,
+        )
+
+    def bipartite_edge_index(
+        self,
+        first_order_graph: Optional[Graph] = None,
+        mapping: str = "last",
+        device: Optional[torch.device] = None,
+    ) -> torch.Tensor:
+        """Return the edge index connecting higher-order nodes to first-order nodes.
+
+        Used by the [DBGNN][pathpyG.nn.dbgnn.DBGNN] model to pass messages from
+        higher-order node representations to first-order ones. Unlike the free function
+        [`generate_bipartite_edge_index`][pathpyG.utils.dbgnn.generate_bipartite_edge_index],
+        this works for any order: "last" refers to the last node of the path, whatever
+        its length.
+
+        Args:
+            first_order_graph: The first-order graph. Optional; accepted so that call
+                sites read symmetrically, and used only for its device.
+            mapping: Which first-order nodes to connect to. One of "last", "first" or "both".
+            device: Device on which to create the tensor.
+
+        Returns:
+            torch.Tensor: Edge index of shape `(2, ·)`, higher-order nodes in the first row.
+        """
+        if device is None:
+            device = first_order_graph.device if first_order_graph is not None else self.device
+
+        node_sequence = self.data.node_sequence
+        ho_idx = torch.arange(self.n, device=device)
+
+        if mapping == "last":
+            fo_idx = node_sequence[:, -1].to(device)
+        elif mapping == "first":
+            fo_idx = node_sequence[:, 0].to(device)
+        elif mapping == "both":
+            fo_idx = torch.cat([node_sequence[:, 0], node_sequence[:, -1]]).to(device)
+            ho_idx = torch.cat([ho_idx, ho_idx])
+        else:
+            raise ValueError(f"Unknown mapping {mapping}. Only 'last', 'first' and 'both' are accepted.")
+
+        return torch.stack([ho_idx, fo_idx])
+
+    @property
+    def n_first_order(self) -> int:
+        """Number of first-order nodes underlying the higher-order nodes."""
+        return self._n_first_order
+
+    def node_id(self, idx: int) -> Union[str, int, tuple]:
+        """Return the first-order path represented by the higher-order node `idx`."""
+        seq = self.data.node_sequence[idx]
+        if self.order == 1:
+            return self.first_order_mapping.to_id(int(seq[0].item()))
+        if self.first_order_mapping.has_ids:
+            return tuple(self.first_order_mapping.to_ids(seq.cpu()).tolist())
+        return tuple(seq.tolist())
+
+    def __str__(self) -> str:
+        """Return a human-readable summary of the higher-order graph."""
+        s = (
+            f"Higher-order graph of order {self.order} with {self.n} nodes and {self.m} edges\n"
+            f"(over {self.n_first_order} first-order nodes)\n"
+        )
+        return s + "\n".join(super().__str__().split("\n")[1:])

From 304fe5779e938c26cf93135b9670cd4670cf7636 Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Thu, 13 Aug 2026 08:06:07 -0400
Subject: [PATCH 03/11] removed some uneeded checks (base constructor already
 does these)

---
 src/pathpyG/core/higher_order_graph.py | 12 ------------
 1 file changed, 12 deletions(-)

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index 56b7bebc..7d3b6a0f 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -6,7 +6,6 @@
 from typing import Optional, Union
 
 import torch
-from torch_geometric import EdgeIndex
 from torch_geometric.data import Data
 from torch_geometric.utils import coalesce
 
@@ -87,14 +86,6 @@ def __init__(
         if "node_sequence" not in data and order not in (None, 1):
             raise ValueError(f"A HigherOrderGraph of order {order} requires a `node_sequence` node attribute.")
 
-        if isinstance(data.edge_index, EdgeIndex):
-            # `Graph.__init__` re-sorts the edge index and reindexes every edge attribute by
-            # the returned permutation - but `EdgeIndex.sort_by` returns `None` for an index
-            # already known to be sorted, and `attr[None]` would add a dimension. Higher-order
-            # graphs are routinely built from already-aggregated (hence sorted) data, so hand
-            # the base class a plain tensor and let it derive a real permutation.
-            data.edge_index = data.edge_index.as_tensor()
-
         super().__init__(data, mapping=mapping)
 
         # `Graph` creates an identity node sequence if none is given, so `self.order`
@@ -195,9 +186,6 @@ def from_aggregated(
         Returns:
             HigherOrderGraph: The aggregated higher-order graph.
         """
-        if isinstance(edge_index, torch.Tensor) and hasattr(edge_index, "as_tensor"):
-            edge_index = edge_index.as_tensor()
-
         order = node_sequence.size(1)
         if first_order_mapping is None:
             first_order_mapping = IndexMap()

From b76281bf6a79cc768da8331e76f3b234c2ef41f8 Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Thu, 13 Aug 2026 08:10:42 -0400
Subject: [PATCH 04/11] removed unused utils.dbgnn module

---
 src/pathpyG/core/higher_order_graph.py |  6 ++--
 src/pathpyG/utils/dbgnn.py             | 46 --------------------------
 tests/nn/test_dbgnn.py                 |  5 ++-
 3 files changed, 4 insertions(+), 53 deletions(-)
 delete mode 100644 src/pathpyG/utils/dbgnn.py

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index 7d3b6a0f..a5301215 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -454,10 +454,8 @@ def bipartite_edge_index(
         """Return the edge index connecting higher-order nodes to first-order nodes.
 
         Used by the [DBGNN][pathpyG.nn.dbgnn.DBGNN] model to pass messages from
-        higher-order node representations to first-order ones. Unlike the free function
-        [`generate_bipartite_edge_index`][pathpyG.utils.dbgnn.generate_bipartite_edge_index],
-        this works for any order: "last" refers to the last node of the path, whatever
-        its length.
+        higher-order node representations to first-order ones. This works for any order:
+        "last" refers to the last node of the path, whatever its length.
 
         Args:
             first_order_graph: The first-order graph. Optional; accepted so that call
diff --git a/src/pathpyG/utils/dbgnn.py b/src/pathpyG/utils/dbgnn.py
deleted file mode 100644
index a70ad8a7..00000000
--- a/src/pathpyG/utils/dbgnn.py
+++ /dev/null
@@ -1,46 +0,0 @@
-"""Utils for DBGNN models."""
-
-from typing import Optional
-
-import torch
-
-from pathpyG.core.graph import Graph
-
-
-def generate_bipartite_edge_index(
-    g: Graph, g2: Graph, mapping: str = "last", device: Optional[torch.device] = None
-) -> torch.Tensor:
-    """Generate edge_index for bipartite graph connecting nodes of a second-order graph to first-order nodes.
-    
-    The mapping strategy determines to which first-order nodes the second-order nodes are connected:
-    - "last": Connects each second-order node to the last node in its sequence.
-    - "first": Connects each second-order node to the first node in its sequence.
-    - "both": Connects each second-order node to both the first and last nodes in its sequence.
-
-    !!! warning "Only for Second-Order Graphs"
-        This function is intended to be used with second-order graphs only. 
-        It does not support the use of higher-order graphs, such as third-order graphs or beyond.
-
-    Args:
-        g (Graph): The first-order graph.
-        g2 (Graph): The second-order graph.
-        mapping (str, optional): The mapping strategy to use. Options are "last", "first", or "both". Defaults to "last".
-        device (torch.device, optional): The device to place the tensor on. Defaults to None.
-
-    Returns:
-        torch.Tensor: The edge_index tensor for the bipartite graph.
-    """
-    if mapping == "last":
-        bipartide_edge_index = torch.tensor([list(range(g2.n)), [v[1] for v in g2.data.node_sequence]], device=device)
-    elif mapping == "first":
-        bipartide_edge_index = torch.tensor([list(range(g2.n)), [v[0] for v in g2.data.node_sequence]], device=device)
-    else:
-        bipartide_edge_index = torch.tensor(
-            [
-                list(range(g2.n)) + list(range(g2.n)),
-                [v[0] for v in g2.data.node_sequence] + [v[1] for v in g2.data.node_sequence],
-            ],
-            device=device,
-        )
-
-    return bipartide_edge_index
diff --git a/tests/nn/test_dbgnn.py b/tests/nn/test_dbgnn.py
index 0de6d8a4..418f8084 100644
--- a/tests/nn/test_dbgnn.py
+++ b/tests/nn/test_dbgnn.py
@@ -5,7 +5,6 @@
 
 from pathpyG.core.multi_order_model import MultiOrderModel
 from pathpyG.nn.dbgnn import DBGNN
-from pathpyG.utils.dbgnn import generate_bipartite_edge_index
 
 
 def test_bipartite_edge_index(simple_walks):
@@ -17,13 +16,13 @@ def test_bipartite_edge_index(simple_walks):
     print(g2.data.edge_index)
     print(g2.mapping)
 
-    bipartite_edge_index = generate_bipartite_edge_index(g, g2, mapping="last")
+    bipartite_edge_index = g2.bipartite_edge_index(g, mapping="last")
     print(bipartite_edge_index)
     # ensure that A,C and B,C are mapped to C, C,D is mapped to D and C,E is mapped to E
 
     assert equal(bipartite_edge_index, tensor([[0, 1, 2, 3], [2, 2, 3, 4]]))
 
-    bipartite_edge_index = generate_bipartite_edge_index(g, g2, mapping="first")
+    bipartite_edge_index = g2.bipartite_edge_index(g, mapping="first")
     print(bipartite_edge_index)
     # ensure that A,C is mapped A, B,C is mapped to B, and C,D and C,E are mapped to C
 

From e2e4f712db863ac5d9b30b6003978461d2d614ef Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Thu, 13 Aug 2026 08:15:01 -0400
Subject: [PATCH 05/11] removed duplicated logic in HigherOrderGraph's order 2
 branch; now delegating to MultiOrderModel

---
 src/pathpyG/core/higher_order_graph.py | 33 +++++++++-----------------
 tests/core/test_event_graph.py         | 30 +++++++++++++++++++++++
 2 files changed, 41 insertions(+), 22 deletions(-)

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index a5301215..85ec9e51 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -338,11 +338,14 @@ def from_path_data(cls, path_data: PathData, order: int = 1, mode: str = "propag
     def from_event_graph(cls, eg: EventGraph, order: int = 2) -> HigherOrderGraph:
         """Aggregate an [`EventGraph`][pathpyG.core.event_graph.EventGraph] into an order-`k` graph.
 
-        For the default order 2, every event whose underlying `(u, v)` pair is the same
-        collapses into a single second-order node, and repeated continuations become an
-        edge weight. Timestamps and the time window `delta` are not represented in the
-        result. Other orders are computed from the time-respecting paths that the event
-        graph encodes, which for orders above 2 means lifting its continuations further.
+        The nodes are the time-respecting paths of `k` first-order nodes that the event
+        graph encodes: events sharing the same underlying path collapse into a single
+        higher-order node, and repeated continuations become an edge weight. Timestamps
+        and the time window `delta` are not represented in the result.
+
+        Equivalent to [`from_temporal_graph`][pathpyG.HigherOrderGraph.from_temporal_graph]
+        on the underlying temporal graph with the event graph's `delta`, but reuses the
+        already-computed continuations instead of lifting the temporal graph again.
 
         Args:
             eg: The second-order temporal event graph to aggregate.
@@ -352,24 +355,10 @@ def from_event_graph(cls, eg: EventGraph, order: int = 2) -> HigherOrderGraph:
             HigherOrderGraph: A higher-order graph of order `order`. It has no nodes if
             there is no time-respecting path of that length.
         """
-        if order != 2:
-            cls._validate_order(order)
-            from pathpyG.core.multi_order_model import MultiOrderModel
-
-            return MultiOrderModel.from_event_graph(eg, max_order=order, cached=False).layers[order]
-
-        edge_index = eg.data.edge_index.as_tensor()
-        # Each continuation carries the weight of the event it starts from, matching the
-        # "src" aggregation used when building order-2 layers from a temporal graph.
-        edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
+        cls._validate_order(order)
+        from pathpyG.core.multi_order_model import MultiOrderModel
 
-        return cls.from_aggregated(
-            edge_index,
-            eg.data.node_sequence,
-            first_order_mapping=eg.first_order_mapping,
-            edge_weight=edge_weight,
-            n_first_order=eg.n_first_order,
-        )
+        return MultiOrderModel.from_event_graph(eg, max_order=order, cached=False).layers[order]
 
     def lift(self, aggr: str = "src") -> HigherOrderGraph:
         """Return the De Bruijn graph of order `k + 1` obtained by lifting this graph.
diff --git a/tests/core/test_event_graph.py b/tests/core/test_event_graph.py
index 94dbf523..039bc734 100644
--- a/tests/core/test_event_graph.py
+++ b/tests/core/test_event_graph.py
@@ -7,6 +7,7 @@
 from torch_geometric.data import Data
 
 from pathpyG.core.event_graph import EventGraph
+from pathpyG.core.higher_order_graph import HigherOrderGraph
 from pathpyG.core.index_map import IndexMap
 from pathpyG.core.multi_order_model import MultiOrderModel
 from pathpyG.core.temporal_graph import TemporalGraph
@@ -218,6 +219,35 @@ def test_multi_order_model_construction(event_graph, temporal_graph):
         )
 
 
+def test_higher_order_graph_from_weighted_event_graph(temporal_graph):
+    """Aggregating an EventGraph respects the edge weights of the temporal graph.
+
+    Regression test: order 2 used to count each continuation once instead of carrying
+    the weight of the event it starts from, so it disagreed with the temporal-graph
+    route for every order but 2.
+    """
+    temporal_graph.data.edge_weight = torch.tensor([2.0, 5.0, 11.0, 7.0])
+    event_graph = EventGraph.from_temporal_graph(temporal_graph, delta=DELTA)
+
+    for k in (1, 2, 3):
+        from_eg = HigherOrderGraph.from_event_graph(event_graph, order=k)
+        from_tg = MultiOrderModel.from_temporal_graph(temporal_graph, delta=DELTA, max_order=k).layers[k]
+
+        assert from_eg.order == k
+        assert from_eg.nodes == from_tg.nodes
+        assert torch.equal(
+            from_eg.data.edge_index.as_tensor(),
+            from_tg.data.edge_index.as_tensor(),
+        )
+        assert torch.equal(from_eg.data.edge_weight, from_tg.data.edge_weight)
+
+    # The weights must actually reflect the temporal graph, not just agree with each other.
+    order_2 = HigherOrderGraph.from_event_graph(event_graph, order=2)
+    assert order_2.nodes == [("a", "b"), ("b", "c"), ("b", "d"), ("c", "e")]
+    # (a,b)->(b,c) carries the weight of event (a->b)@1, (b,c)->(c,e) that of (b->c)@2
+    assert order_2.data.edge_weight.tolist() == [2.0, 5.0]
+
+
 def test_to_device(event_graph):
     """Moving an EventGraph moves its underlying TemporalGraph too."""
     moved = event_graph.to(torch.device("cpu"))

From 018b6831869ad18bbb651cd0e56e265cba11991f Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Thu, 13 Aug 2026 08:18:26 -0400
Subject: [PATCH 06/11] added two helper functions to algorithms to reduce code
 duplication

---
 src/pathpyG/algorithms/lift_order.py   | 51 ++++++++++++++++++++++++++
 src/pathpyG/core/higher_order_graph.py |  9 ++---
 src/pathpyG/core/multi_order_model.py  | 15 +++-----
 3 files changed, 60 insertions(+), 15 deletions(-)

diff --git a/src/pathpyG/algorithms/lift_order.py b/src/pathpyG/algorithms/lift_order.py
index 1b2269b9..1d95e4b2 100644
--- a/src/pathpyG/algorithms/lift_order.py
+++ b/src/pathpyG/algorithms/lift_order.py
@@ -106,6 +106,57 @@ def lift_order_edge_index_weighted(
     return ho_index, ho_edge_weight
 
 
+def lift_node_sequence(edge_index: torch.Tensor, node_sequence: torch.Tensor) -> torch.Tensor:
+    """Extend node sequences by one order along an edge index.
+
+    Each edge `(u, v)` of the (k-1)-th order graph becomes a node of the k-th order graph,
+    representing the path of `u` followed by the last first-order node of `v`.
+
+    Args:
+        edge_index: A **sorted** edge index tensor of shape (2, num_edges).
+        node_sequence: The node sequences of the (k-1)-th order graph, of shape (num_nodes, k-1).
+
+    Returns:
+        The node sequences of the k-th order graph, of shape (num_edges, k).
+    """
+    return torch.cat([node_sequence[edge_index[0]], node_sequence[edge_index[1]][:, -1:]], dim=1)
+
+
+def lift_order_step(
+    edge_index: torch.Tensor,
+    node_sequence: torch.Tensor,
+    edge_weight: torch.Tensor | None = None,
+    aggr: str = "src",
+) -> tuple[torch.Tensor, torch.Tensor, torch.Tensor | None]:
+    """Lift an edge index together with its node sequences by one order.
+
+    Combines the line-graph transformation of the edge index with the corresponding
+    extension of the node sequences, so that the result again describes a graph whose
+    nodes are paths of first-order nodes. The result is **not** aggregated: duplicate
+    node sequences are left for [`aggregate_edge_index`][pathpyG.algorithms.lift_order.aggregate_edge_index]
+    (or [`HigherOrderGraph.from_aggregated`][pathpyG.HigherOrderGraph.from_aggregated]) to merge.
+
+    Args:
+        edge_index: A **sorted** edge index tensor of shape (2, num_edges).
+        node_sequence: The node sequences of the (k-1)-th order graph.
+        edge_weight: The edge weights of the (k-1)-th order graph. If None, the lifted
+            graph is returned without weights.
+        aggr: The aggregation method for the edge weights. One of "src", "dst", "max",
+            "mul" or "add". Ignored if `edge_weight` is None.
+
+    Returns:
+        A tuple of the lifted edge index, the lifted node sequences and the aggregated
+        edge weights (None if `edge_weight` was None).
+    """
+    if edge_weight is None:
+        ho_index = lift_order_edge_index(edge_index, num_nodes=node_sequence.size(0))
+    else:
+        ho_index, edge_weight = lift_order_edge_index_weighted(
+            edge_index, edge_weight=edge_weight, num_nodes=node_sequence.size(0), aggr=aggr
+        )
+    return ho_index, lift_node_sequence(edge_index, node_sequence), edge_weight
+
+
 def aggregate_edge_index(
     edge_index: torch.Tensor, node_sequence: torch.Tensor, edge_weight: torch.Tensor | None = None, aggr: str = "sum"
 ) -> Graph:
diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index 85ec9e51..1403bd5a 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -9,7 +9,7 @@
 from torch_geometric.data import Data
 from torch_geometric.utils import coalesce
 
-from pathpyG.algorithms.lift_order import aggregate_edge_index, lift_order_edge_index_weighted
+from pathpyG.algorithms.lift_order import aggregate_edge_index, lift_order_step
 from pathpyG.core.event_graph import EventGraph
 from pathpyG.core.graph import Graph
 from pathpyG.core.index_map import IndexMap
@@ -385,11 +385,8 @@ def lift(self, aggr: str = "src") -> HigherOrderGraph:
         else:
             edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
 
-        ho_index, ho_weight = lift_order_edge_index_weighted(
-            edge_index, edge_weight=edge_weight, num_nodes=self.n, aggr=aggr
-        )
-        node_sequence = torch.cat(
-            [self.data.node_sequence[edge_index[0]], self.data.node_sequence[edge_index[1]][:, -1:]], dim=1
+        ho_index, node_sequence, ho_weight = lift_order_step(
+            edge_index, self.data.node_sequence, edge_weight=edge_weight, aggr=aggr
         )
 
         return HigherOrderGraph.from_aggregated(
diff --git a/src/pathpyG/core/multi_order_model.py b/src/pathpyG/core/multi_order_model.py
index 88d068f7..1da9fd33 100644
--- a/src/pathpyG/core/multi_order_model.py
+++ b/src/pathpyG/core/multi_order_model.py
@@ -13,8 +13,9 @@
 from pathpyG.algorithms.lift_order import (
     aggregate_edge_index,
     aggregate_node_attributes,
+    lift_node_sequence,
     lift_order_edge_index,
-    lift_order_edge_index_weighted,
+    lift_order_step,
 )
 from pathpyG.core.event_graph import EventGraph
 from pathpyG.core.higher_order_graph import HigherOrderGraph
@@ -115,13 +116,9 @@ def iterate_lift_order(
                 Defaults to the number of IDs in `mapping`.
         """
         # Lift order
-        if edge_weight is None:
-            ho_index = lift_order_edge_index(edge_index, num_nodes=node_sequence.size(0))
-        else:
-            ho_index, edge_weight = lift_order_edge_index_weighted(
-                edge_index, edge_weight=edge_weight, num_nodes=node_sequence.size(0), aggr=aggr
-            )
-        node_sequence = torch.cat([node_sequence[edge_index[0]], node_sequence[edge_index[1]][:, -1:]], dim=1)
+        ho_index, node_sequence, edge_weight = lift_order_step(
+            edge_index, node_sequence, edge_weight=edge_weight, aggr=aggr
+        )
 
         # Aggregate
         if save:
@@ -180,7 +177,7 @@ def from_temporal_graph(
             )
 
         if max_order > 1:
-            node_sequence = torch.cat([node_sequence[edge_index[0]], node_sequence[edge_index[1]][:, -1:]], dim=1)
+            node_sequence = lift_node_sequence(edge_index, node_sequence)
             if event_graph is None:
                 edge_index = EventGraph.build_edge_index(g, delta)
             else:

From 2768cfea808083161f04d8f05715a0921fbdd8bf Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Thu, 13 Aug 2026 09:51:17 -0400
Subject: [PATCH 07/11] removed unhelpful comment

---
 src/pathpyG/core/higher_order_graph.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index 1403bd5a..7fb7fc08 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -302,7 +302,6 @@ def from_temporal_graph(
             [('a', 'c'), ('c', 'd')]
         """
         cls._validate_order(order)
-        # Imported here because `MultiOrderModel` builds `HigherOrderGraph` layers itself.
         from pathpyG.core.multi_order_model import MultiOrderModel
 
         return MultiOrderModel.from_temporal_graph(

From cb7bd2febc02e0c00b4c6818a148a8392a533e42 Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Thu, 13 Aug 2026 10:07:43 -0400
Subject: [PATCH 08/11] comments

---
 src/pathpyG/core/higher_order_graph.py | 34 +++-----------------------
 src/pathpyG/core/multi_order_model.py  |  3 +--
 2 files changed, 5 insertions(+), 32 deletions(-)

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index 7fb7fc08..2cb79773 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -22,11 +22,9 @@
 class HigherOrderGraph(Graph):
     """A De Bruijn graph of order `k`, whose nodes are paths of `k` first-order nodes.
 
-    Where a [`Graph`][pathpyG.Graph] has one node per entity and an
-    [`EventGraph`][pathpyG.core.event_graph.EventGraph] has one node per observed
-    interaction, a `HigherOrderGraph` has one node per *distinct* path of length `k`
+    A `HigherOrderGraph` has one node per distinct path of length `k`
     in the underlying first-order graph. Repeated observations of the same path are
-    aggregated into an `edge_weight`, so timestamps are no longer represented: this
+    aggregated into an `edge_weight`. Timestamps are not represented: this
     is a model of how paths flow rather than a record of what happened.
 
     Order 1 is the degenerate case and is simply the weighted first-order graph, with
@@ -170,10 +168,6 @@ def from_aggregated(
     ) -> HigherOrderGraph:
         """Aggregate a (possibly duplicated) higher-order edge index into a De Bruijn graph.
 
-        This is the single place where higher-order nodes get their identity: duplicate
-        node sequences are merged, edge weights are aggregated, and the higher-order
-        `IndexMap` naming each node by its path is built.
-
         Args:
             edge_index: Edge index whose nodes are indices into `node_sequence`.
             node_sequence: Tensor of shape `(num_nodes, order)` with the first-order path
@@ -220,9 +214,6 @@ def from_aggregated_graph(
     ) -> HigherOrderGraph:
         """Adopt an already-aggregated [`Graph`][pathpyG.Graph] as a higher-order graph.
 
-        Used to give the layers computed by a multi-order model their proper type. The
-        underlying `data` object is shared, not copied.
-
         Args:
             g: Aggregated graph carrying a `node_sequence` of shape `(num_nodes, order)`.
             first_order_mapping: Mapping of the underlying first-order node IDs.
@@ -282,8 +273,7 @@ def from_temporal_graph(
         Order 1 is simply the weighted static graph and ignores `delta`; for higher orders
         the nodes are the time-respecting paths of `k` nodes, i.e. those whose consecutive
         interactions are at most `delta` apart. Orders above 2 are reached by repeatedly
-        lifting the *unaggregated* data, so the edge weights count observed paths rather
-        than being implied by lower-order statistics (unlike [`lift`][pathpyG.HigherOrderGraph.lift]).
+        lifting the unaggregated data.
 
         Args:
             g: The temporal graph.
@@ -337,14 +327,8 @@ def from_path_data(cls, path_data: PathData, order: int = 1, mode: str = "propag
     def from_event_graph(cls, eg: EventGraph, order: int = 2) -> HigherOrderGraph:
         """Aggregate an [`EventGraph`][pathpyG.core.event_graph.EventGraph] into an order-`k` graph.
 
-        The nodes are the time-respecting paths of `k` first-order nodes that the event
-        graph encodes: events sharing the same underlying path collapse into a single
-        higher-order node, and repeated continuations become an edge weight. Timestamps
-        and the time window `delta` are not represented in the result.
-
         Equivalent to [`from_temporal_graph`][pathpyG.HigherOrderGraph.from_temporal_graph]
-        on the underlying temporal graph with the event graph's `delta`, but reuses the
-        already-computed continuations instead of lifting the temporal graph again.
+        on the underlying temporal graph with the event graph's `delta`.
 
         Args:
             eg: The second-order temporal event graph to aggregate.
@@ -365,12 +349,6 @@ def lift(self, aggr: str = "src") -> HigherOrderGraph:
         Nodes of the result are the edges of this graph, i.e. the paths of length `k + 1`
         that exist in this graph's topology.
 
-        Warning:
-            This lifts an *aggregated* graph, so the resulting edge weights are those
-            implied by the order-`k` statistics rather than counts of observed paths of
-            length `k + 1`. To fit a layer to observations, use
-            [`MultiOrderModel`][pathpyG.MultiOrderModel], which lifts the unaggregated data.
-
         Args:
             aggr: Aggregation used for the lifted edge weights. One of "src", "dst",
                 "max", "mul" or "add".
@@ -438,10 +416,6 @@ def bipartite_edge_index(
     ) -> torch.Tensor:
         """Return the edge index connecting higher-order nodes to first-order nodes.
 
-        Used by the [DBGNN][pathpyG.nn.dbgnn.DBGNN] model to pass messages from
-        higher-order node representations to first-order ones. This works for any order:
-        "last" refers to the last node of the path, whatever its length.
-
         Args:
             first_order_graph: The first-order graph. Optional; accepted so that call
                 sites read symmetrically, and used only for its device.
diff --git a/src/pathpyG/core/multi_order_model.py b/src/pathpyG/core/multi_order_model.py
index 1da9fd33..6f75f843 100644
--- a/src/pathpyG/core/multi_order_model.py
+++ b/src/pathpyG/core/multi_order_model.py
@@ -34,8 +34,7 @@ class MultiOrderModel:
     Each graph layer is represented as a
     [HigherOrderGraph][pathpyG.core.higher_order_graph.HigherOrderGraph] object, layer 1
     included. Each layer therefore knows its own order and the first-order nodes it was
-    built from, so results can be projected back onto entities without the caller keeping
-    track of it.
+    built from.
     This class provides methods to search for the optimal order of the model based on likelihood ratio tests,
     as well as methods to compute the log-likelihood of observed paths given the model.
 

From 4112444d0ca6b0620ca2923c1e5f89d860d5bb07 Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Thu, 13 Aug 2026 10:30:03 -0400
Subject: [PATCH 09/11] explanatory note

---
 src/pathpyG/core/higher_order_graph.py | 15 +++++++++++++++
 1 file changed, 15 insertions(+)

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index 2cb79773..b856a0c2 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -285,6 +285,11 @@ def from_temporal_graph(
             HigherOrderGraph: A higher-order graph of order `order`. It has no nodes if
             there is no time-respecting path of that length.
 
+        Note:
+            Each call rebuilds the whole chain of lifts from order 1. To obtain several
+            orders, build a [`MultiOrderModel`][pathpyG.MultiOrderModel] with
+            `cached=True` once and read its `layers` instead.
+
         Examples:
             >>> import pathpyG as pp
             >>> t = pp.TemporalGraph.from_edge_list([("a", "c", 1), ("c", "d", 2)])
@@ -311,6 +316,11 @@ def from_path_data(cls, path_data: PathData, order: int = 1, mode: str = "propag
             HigherOrderGraph: A higher-order graph of order `order`. It has no nodes if
             no observed path is that long.
 
+        Note:
+            Each call rebuilds the whole chain of lifts from order 1. To obtain several
+            orders, build a [`MultiOrderModel`][pathpyG.MultiOrderModel] with
+            `cached=True` once and read its `layers` instead.
+
         Examples:
             >>> import pathpyG as pp
             >>> paths = pp.PathData(pp.IndexMap(list("acd")))
@@ -337,6 +347,11 @@ def from_event_graph(cls, eg: EventGraph, order: int = 2) -> HigherOrderGraph:
         Returns:
             HigherOrderGraph: A higher-order graph of order `order`. It has no nodes if
             there is no time-respecting path of that length.
+
+        Note:
+            Each call rebuilds the whole chain of lifts from order 1. To obtain several
+            orders, build a [`MultiOrderModel`][pathpyG.MultiOrderModel] with
+            `cached=True` once and read its `layers` instead.
         """
         cls._validate_order(order)
         from pathpyG.core.multi_order_model import MultiOrderModel

From 809b49581a0935088573fe7a1210dbf5b21d774a Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Wed, 2 Sep 2026 12:17:18 -0400
Subject: [PATCH 10/11] from_aggregated -> aggregate renaming

---
 src/pathpyG/core/higher_order_graph.py | 6 +++---
 src/pathpyG/core/multi_order_model.py  | 6 +++---
 2 files changed, 6 insertions(+), 6 deletions(-)

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index b856a0c2..997dcbe1 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -157,7 +157,7 @@ def _build_mapping(node_sequence: torch.Tensor, first_order_mapping: IndexMap) -
         return IndexMap([tuple(v.tolist()) for v in node_sequence])
 
     @classmethod
-    def from_aggregated(
+    def aggregate(
         cls,
         edge_index: torch.Tensor,
         node_sequence: torch.Tensor,
@@ -252,7 +252,7 @@ def from_graph(cls, g: Graph, weight: str = "edge_weight") -> HigherOrderGraph:
             edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
         node_sequence = torch.arange(g.n, device=edge_index.device).unsqueeze(1)
 
-        return cls.from_aggregated(
+        return cls.aggregate(
             edge_index,
             node_sequence,
             first_order_mapping=g.mapping,
@@ -381,7 +381,7 @@ def lift(self, aggr: str = "src") -> HigherOrderGraph:
             edge_index, self.data.node_sequence, edge_weight=edge_weight, aggr=aggr
         )
 
-        return HigherOrderGraph.from_aggregated(
+        return HigherOrderGraph.aggregate(
             ho_index,
             node_sequence,
             first_order_mapping=self.first_order_mapping,
diff --git a/src/pathpyG/core/multi_order_model.py b/src/pathpyG/core/multi_order_model.py
index 6f75f843..30d02a1d 100644
--- a/src/pathpyG/core/multi_order_model.py
+++ b/src/pathpyG/core/multi_order_model.py
@@ -121,7 +121,7 @@ def iterate_lift_order(
 
         # Aggregate
         if save:
-            gk = HigherOrderGraph.from_aggregated(
+            gk = HigherOrderGraph.aggregate(
                 ho_index,
                 node_sequence,
                 first_order_mapping=mapping,
@@ -167,7 +167,7 @@ def from_temporal_graph(
         else:
             edge_weight = torch.ones(edge_index.size(1), device=edge_index.device)
         if cached or max_order == 1:
-            m.layers[1] = HigherOrderGraph.from_aggregated(
+            m.layers[1] = HigherOrderGraph.aggregate(
                 edge_index=edge_index,
                 node_sequence=node_sequence,
                 edge_weight=edge_weight,
@@ -185,7 +185,7 @@ def from_temporal_graph(
 
             # Aggregate
             if cached or max_order == 2:
-                m.layers[2] = HigherOrderGraph.from_aggregated(
+                m.layers[2] = HigherOrderGraph.aggregate(
                     edge_index=edge_index,
                     node_sequence=node_sequence,
                     edge_weight=edge_weight,

From 8d9bd34d03ffd0135962254165ac9d429037935b Mon Sep 17 00:00:00 2001
From: Vineet Bansal <vineetbansal@protonmail.com>
Date: Tue, 8 Sep 2026 12:51:56 -0400
Subject: [PATCH 11/11] Update src/pathpyG/core/higher_order_graph.py

Co-authored-by: Moritz Lampert <moritz@lampert.ai>
---
 src/pathpyG/core/higher_order_graph.py | 9 +++++++++
 1 file changed, 9 insertions(+)

diff --git a/src/pathpyG/core/higher_order_graph.py b/src/pathpyG/core/higher_order_graph.py
index 997dcbe1..378f708d 100644
--- a/src/pathpyG/core/higher_order_graph.py
+++ b/src/pathpyG/core/higher_order_graph.py
@@ -395,6 +395,15 @@ def to_first_order(self, mode: str = "last") -> Graph:
         Each higher-order node is replaced by one of the first-order nodes of its path,
         and the weights of higher-order edges mapping to the same first-order edge are
         summed. First-order nodes not traversed by any path remain as isolated nodes.
+        
+        Warning: This is a projection, not an inverse transformation
+            This method does not reconstruct the original first-order graph from
+            which this higher-order graph was built. Instead, it maps each higher-
+            order node to either the first or last first-order node in its represented path.
+            
+            Consequently, the result preserves flow encoded by the higher-order model
+            under the selected projection, but may differ from the original graph in its
+            edge set. In particular, isolated first-order edges cannot be recovered.
 
         Args:
             mode: Which first-order node of the path represents it. Either "last" or "first".