PCEvo: Path-Consistent Molecular Representation via Virtual Evolutionary

1 Introduction

Recent years have witnessed rapid progress in molecular representation learning driven by LLMs-based [Zheng et al., 2025; Yu et al., 2025b; Li et al., 2025a], graph-based [Yu et al., 2025a; Yu et al., 2024], and geometry-aware neural architectures [Rao et al., 2026; Lin et al., 2026; Xiong et al., 2025]. Message passing networks and 3D geometric models have become mainstream for learning structuredependent signals [Qin et al., 2026; Li et al., 2025b], with representative backbones spanning MPNNs [Gilmer et al., 2017], SchNet [Sch¨utt et al., 2018], DimeNet [Gasteiger et al., 2020], and more recent equivariant Transformer families such as Equiformer [Liao and Smidt, 2023] and ViSNet [Wang et al., 2024b]. In parallel, self-supervised and contrastive pretraining has been increasingly adopted to leverage large unlabeled corpora and improve data efficiency [Li et al., 2025c], including general graph contrastive objectives [Wang et al., 2024a] and molecule-specific pretraining paradigms such as GROVER [Rong et al., 2020], MolCLR [Wang et al., 2022], and UniMol [Zhou et al., 2023]. These developments are routinely evaluated on established benchmarks such as QM9 [Ramakrishnan et al., 2014] and MoleculeNet [Wu et al., 2018], and are further expanding toward crossmodal alignment between molecular structures and natural language descriptions [Wang et al., 2025].

Despite this progress, the dominant learning paradigm remains largely static: each molecule is treated as a single endpoint sample (i.e., only its final observed structure is used), and the model is trained to regress properties via direct supervision on sparse labels [An et al., 2025; Zhang et al., 2025] (as shown in Fig. 1). In few-shot settings, this formulation is prone to overfitting because chemical space is sparsely covered; as a result, models can rely on dataset-specific correlations and frequent motifs, leading to unstable generalization. More fundamentally, most static molecular representations do not explicitly model how properties change under chemically feasible local edits to molecular structure [Li et al., 2026]. In practical drug discovery, where labeled assays are costly and task distributions are sparse, capturing such regularities can improve the robustness of learned representations under the few-shot setting.

From a chemical perspective, molecules are discrete combinatorial objects whose structures can be naturally described as graphs and encoded by rule-based line notations [Gilmer et al., 2017; Weininger, 1988]. More importantly, property changes in medicinal chemistry are often analyzed through small, chemically valid transformations, where a minimal structural modification induces a measurable shift in activity or physicochemical profiles [Ishikura et al., 2025; Hussain and Rea, 2010a]. These observations motivate a virtual evolutionary path view of the structure-property relationship (Fig. 2) [Maggiora and Bajorath, 2014; Stumpfe and Bajorath, 2012]. By treating minimal edit operations, such as atom substitution, local connectivity adjustment, and functional group modification, as fundamental units [Walz et al., 2025; Hussain and Rea, 2010b], we can construct virtual evolutionary paths connecting different molecular states [Jensen, 2019]. Modeling property variation along such paths enables the models to leverage incremental structural signals that are often underutilized when training solely on static endpoint samples.

Based on this perspective, we propose PCEvo, a pathconsistent representation method via virtual evolutionary. PCEvo does not treat molecules as isolated, static samples; instead, it represents the structural differences between two similar molecules as an evolutionary sequence of chemically feasible minimal edit operations. PCEvo further converts the labels of the two endpoint molecules into stepwise supervision distributed along the path, enabling the model to learn how incremental structural changes drive property variations during training, rather than performing the regression only at the endpoints. Importantly, there are typically multiple feasible editing paths connecting the same molecular pair; yet, the cumulative property change should remain consistent regardless of which path is taken. Leveraging this principle, PCEvo introduces a path-consistency objective that aligns the accumulated predictions across alternative feasible paths. This objective integrates seamlessly with standard static supervised training, thereby improving generalization in few-shot settings. We evaluate PCEvo on QM9 and MoleculeNet by integrating it into a range of representative backbone methods, thereby systematically demonstrating its generality and effectiveness; results show that PCEvo consistently reduces prediction error under few-shot settings on QM9 and achieves state-of-the-art (SOTA) performance on the three MoleculeNet regression tasks under the standard split, while maintaining reliable gains in both accuracy and stability when labeled data are limited. Our contributions are summarized as follows:

• We propose PCEvo, a path-consistent molecular representation method that learns from virtual evolutionary edit paths between similar molecules, transforming sparse end-point supervision into stepwise training signals along chemically feasible paths.
• Extensive experiments on QM9 and MoleculeNet under few-shot settings show that PCEvo consistently reduces prediction error and improves predictive stability, achieving SOTA performance compared with strong baselines.

2 Methods

2.1 Overview

To address the generalization limitations of static end-to-end learning in few-shot scenarios, as illustrated in Fig. 3, we propose a path-consistent molecular representation method, PCEvo. Deviating from conventional static structural mapping paradigms, PCEvo decomposes complex structureproperty correlations into a sequence of learnable, cumulative edit steps. By imposing these path constraints, PCEvo could capture the intrinsic logic underlying molecular evolution.

2.2 Virtual Evolutionary Path Generation

To transform static molecular graphs into dynamic evolutionary paths, we propose a virtual evolutionary path generation pipeline. This pipeline converts a standard dataset D = {(Gi, yi)}Ni=1 into an evolutionary dataset Devo = {(Gs,Gt, τ) | τ ∈ T (Gs,Gt)}. For each molecular pair (Gs,Gt), T (Gs,Gt) denotes the set of all valid evolutionary paths connecting the source molecule Gs to the target molecule Gt. Each path τ ∈ T (Gs,Gt) is a sequence of edit operations. Explicitly modeling multiple valid paths for each molecular pair allows the PCEvo to account for the intrinsic diversity of structural transformations, facilitating the

1 - Virtual Evolutionary Path Generation

2 - Differential Encoding Unit

learning of property changes that are invariant to the ordering of operations.

Structural Neighbor Retrieval

Constructing evolutionary paths between arbitrary pairs of molecules is chemically implausible and computationally inefficient. To model realistic structural optimizations, we restrict the source molecules to the chemical neighborhood of the target.

For each target molecule Gt in the training or validation set, we query a candidate pool Gs to identify the topK nearest neighbors. We utilize the Tanimoto similarity fTanimoto(Gs,Gt) coefficient [Bajusz et al., 2015] based on extended-connectivity fingerprints as the metric:

This retrieval-based pairing strategy ensures that Gs and Gt share a significant structural scaffold, thereby defining a chemically meaningful optimization landscape.

Fundamental Path Operation Edit Units

Identifying the minimal set of edit operations ensures that the molecular transformation includes only essential structural changes, improving path generation efficiency and facilitating the learning of structure-property relationships. Given a pair (Gs,Gt), we seek the minimal set of edit operations S required to transform Gs into Gt. We formulate this as a maximum common subgraph (MCS) problem [Fuchs and Riesen, 2025; Raymond and Willett, 2002; Bunke, 1997].

We establish an atom-level alignment π : Vs → Vt by minimizing the matching cost matrix defined by atomic numbers and valence constraints. Based on the atom-level alignment π, the symmetric difference between the two graphs is decomposed into a set of atomic operations. The definitions of these operations are detailed in Table 1. Specifically, unmapped atoms in the source molecule Gs correspond to REMOVE ATOM operations, while unmapped atoms in the target molecule Gt correspond to ADD ATOM operations. At the bond level, discrepancies between mapped atom pairs (u, v) in Gs and (π(u), π(v)) in Gt give rise to ADD BOND, REMOVE BOND, or CHANGE BOND operations.

Applying this procedure yields an unordered set of operations S = {o1, . . . , oM}, which constitutes the minimal set of edits required to transformGs into Gt.

Topologically Constrained Path Sampling

While the set S specifies the structural changes, it does not determine the order in which they should be applied. Performing the operations in an arbitrary sequence may violate topological constraints, such as adding a bond to an atom that does not yet exist. To resolve this, we construct a directed acyclic graph, termed the dependency graph: Gdep = (S, Edep). An edge oi → oj exists if operation oi is a topological prerequisite for oj. The core constraints include:

• Existence precondition: An ADD BOND(u, v) operation depends on the ADD ATOM operations that create nodes u and v.
• Deletion latency: A REMOVE ATOM(u) operation must succeed all bond removals or modifications involving node u.

Any linear extension (topological sort) of Gdep constitutes a valid path τ. Importantly, a single edit set S often allows

Algorithm 1 Virtual Evolutionary Path Generation

Require: Standard dataset D, Number of neighbors K, Max

paths Pmax Ensure: Evolutionary dataset Devo 1: Initialize Devo ← ∅ 2: for all Gt ∈ D do 3: Candidate pool Gs ← (Gt, Gs,K) 4: for all Gs ∈ Gs do 5: π ← ALIGNMENT(Gs, Gt) 6: S ← EXTRACTMINIMALEDITS(Gs, Gt, π) 7: Gdep ← BUILDDEPENDENCYGRAPH(S) 8: Tpair ← ∅ 9: while |Tpair| < Pmax do 10: τ ← SAMPLING(Gdep) 11: Tpair ← Tpair ∪ {τ } 12: end while 13: for all τ ∈ Tpair do 14: Devo ← Devo ∪ {(Gs, Gt, τ )} 15: end for 16: end for 17: end for 18: return Devo

exponentially many valid permutations. We leverage this by randomly sampling up to Pmax distinct topological sorts for each molecular pair. This serves as a powerful combinatorial data augmentation strategy: the model learns to associate the same net property change ∆y with multiple distinct paths {τ1, . . . , τPmax }, thereby enforcing path-independence of the physical property and improving generalization. Each path τ ∈ T (Gs, Gt) is a sequence of edit operations, and its length may vary across different paths.

2.3 Differential Evolutionary Path Encoder

To represent the feature changes induced by each virtual edit operation along a molecular structural evolution path, a differential evolutionary path encoder is proposed. Unlike the standard paradigm that encodes molecules as isolated static data points, this module explicitly models the incremental changes between consecutive intermediate states along the virtual path τ .

Let G(t) s denote the intermediate molecular state after the t-th edit operation in the path, with the boundary conditions G(0) s = Gs and G(T) s = Gt, where T is the total path length.

To bridge the discrete symbolic space of chemical edits and the continuous latent space, we first define a bijective mapping ψ : O → RD . The t-th edit operation ot, which triggers the transition G(t) s → G(t+1) s , is vectorized via concatenated encodings:

where components encode the operation type, atomic semantics, bond attributes, and positional indices, respectively. This vectorization provides the initial condition for the specific structural transformation at step t.

Then, we employ a shared-weight molecular encoder fθ(·) to project the adjacent states G(t) s and G(t+1) s into a continuous latent manifold:

Our method PCEvo is architecture-agnostic and shared across all steps to ensure consistent feature extraction along the path. For each virtual edit, we model its incremental effect by differencing adjacent representations and predicting the corresponding property change:

which learns edit-wise contributions in an additive manner.

2.4 Learning Objective

Our optimization jointly anchors absolute property prediction in the global chemical space and enforces summation consistency of property variations along virtual evolutionary paths. The model consists of a shared-weight backbone encoder fθ(·) and a non-linear delta predictor ϕψ(·). Given a labeled molecular pair (Gs,Gt) with properties (ys, yt) and a feasible edit path τ = (G(0) s , . . . ,G(T) s ) satisfyingG(0) s = Gs and G(T ) s = Gt, we denote the latent representation at step t as

Static Property Loss. We supervise absolute predictions at the observed endpoints:

Path Summation Consistency Loss. We require that the cumulative sum of predicted stepwise property increments along the path matches the ground-truth endpoint difference:

where the non-linear predictor ϕψ(·) ensures the summation does not mathematically collapse into a simple difference of endpoints. This forces the gradient signal to propagate through every intermediate state, compelling the model to learn the specific physical contribution of each elementary edit operation. The final objective is L = Lstatic + Lcons.

3 Theoretical Analysis

This section explains why PCEvo generalizes better by decomposing sparse static endpoint regression supervision into edit-level signals along evolutionary paths. Two mechanisms are emphasized: (i) edit decomposition increases the effective number of supervised signals; (ii) learning on edit-local representation increments reduces hypothesis complexity.

Table 4: Overall experiments on MoleculeNet under different train/test splits. RMSE is reported for each task, where R denotes the training sample ratio.

3.1 Learning Formulation

Given a labeled dataset D = {(Gi, yi)}Ni=1 and the induced evolutionary dataset Devo, we compare static endpoint learning and path-consistent learning. Here ̂R(·) denotes the empirical risk (sample average of the loss) minimized during training.

Static model. The static objective Lstatic supervises absolute property prediction on D. We denote by Hstatic the hypothesis class of all static predictors realizable by the backbone when parameters vary:

Path-consistent learning. The consistency objective Lcons enforces that the summed stepwise deltas along τ = (G(0) s , . . . ,G(T) s ) match the endpoint difference (yt - ys). With ht = fθ(G(t) s ), the predicted total change is

We denote by Hcons the hypothesis class of all delta predictors realizable by the MLP head when its parameters vary:

3.2 Generalization Bound

We adopt Rademacher complexity to compare the generalization of static endpoint learning and edit-path learning.

Assumption 1 (Effective sample size under mixing). Editlevel samples induced by virtual paths form a weakly dependent sequence. For instance, they satisfy β-mixing with

∑k≥1 β(k) < ∞. Then the edit-level learning admits an effective sample size

where ¯L is the average path length.

Assumption 1 states that, although edit-level samples generated along paths are weakly dependent, their statistical effect can be captured by an effective sample size neff, so the generalization bound applies by replacing n with neff .

Let Rn(H) denote the empirical Rademacher complexity. By the standard Rademacher generalization theory [Yin et al., 2019; Bartlett et al., 2005], the generalization gap scales as O(Rn(H) + √log(1/δ)/n) uniformly over f ∈ H. Under Assumption 1, edit-level learning obeys the same scaling with n replaced by neff (up to constants absorbed into cmix).

Theorem 1 (Edit decomposition yields a tighter generalization bound). Assume ℓ(·, ·) is bounded in [0, 1] and is 1Lipschitz in its first argument. Assume ∥h∥2 ≤ Bstatic and ∥∆h∥2 ≤ Bedit. If hypotheses in Hstatic and Hcons are LLipschitz with respect to their inputs, then:

(i) Complexity scaling. The empirical Rademacher complexity scales with the input radius and the sample size. In particular,

(ii) Generalization bounds. With probability at least 1-δ,

The two bounds share the same form, but differ in the effective sample size and the input radius. Edit-path learning replaces (N,Bstatic) with (neff,Bedit), where neff is induced by path decomposition and ∆h is edit-local. In typical settings, a single edit produces a local perturbation in representation space, so Bedit ≤ Bstatic. Therefore, the edit-path bound becomes tighter whenever

In particular, the static generalization gap scales as

while edit-path learning follows the same scaling with (N,Bstatic) replaced by (neff,Bedit). Finally, enforcing path consistency further restricts the admissible solutions within Hcons, which can only reduce the effective complexity and tighten the above bound.

4 Experiments

Datasets. We conduct our evaluations on the QM9 dataset [Ramakrishnan et al., 2014], focusing on the prediction of the highest occupied molecular orbital (ϵHOMO) and lowest unoccupied molecular orbital (ϵLUMO). For physicochemical property prediction, we further employ the ESOL, FreeSolv, and Lipophilicity datasets fromMoleculeNet [Wu et al., 2018].

Baselines. We evaluate PCEvo against a comprehensive set of backbone methods, grouped into three categories: 2D GNNs/Transformer: GCN [Kipf and Welling, 2017], GraphSAGE [Hamilton et al., 2017], GAT [Veliˇckovi´c et al., 2018], GIN [Xu et al., 2019], and Moleculeformer [Qin et al., 2025]; 3D geometry-aware: SchNet [Sch¨utt et al., 2018], DimeNet [Gasteiger et al., 2020], Equiformer [Liao and Smidt, 2023], ViSNet [Wang et al., 2024b], and GotenNet [Aykent and Xia, 2025]; Pretraining: N-GRAM [Liu et al., 2019], GROVER [Rong et al., 2020], MolCLR [Wang et al., 2022], Hu et al. [Hu et al., 2020], and EMPP [An et al., 2025].

Implementation Details. We designed specific data splitting protocols to simulate varying degrees of data scarcity. For the QM9 dataset, we constructed extremely few-shot training/testing sets containing randomly sampled 100 and 1,000 molecules. For the smaller MoleculeNet datasets, we simulated few-shot learning by using 10%, 20%, and 80% of the total data for training; the 10% data were allocated for testing. Performance was assessed using three standard regression metrics: mean absolute error (MAE), mean squared error (MSE), and the Pearson correlation coefficient (PCC).

4.1 Overall Experiments

To evaluate the effectiveness and scalability of PCEvo for molecular representation, we conduct experiments on QM9 and MoleculeNet. Following our path construction procedure, for each target molecule, we retrieve the top-N=5 structural neighbors and enumerate up to Pmax=50 valid edit paths for each molecular pair under topological dependency constraints.

On QM9, we adopt two low-resource settings with N=100 and N=1000 labeled training samples (results in Table 2). PCEvo yields consistent improvements across four 3D backbones (SchNet, DimeNet, Equiformer, and ViSNet), and the correlation recovery is particularly pronounced in few-shot settings (N=100): for ϵHOMO, Equiformer improves PCC from 0.07 to 0.53 (+0.46). Moreover, compared with recent strong baselines such as EMPP, Moleculeformer, and GotenNet, PCEvo achieves better performance under the same fewshot protocol. These results indicate that when the chemical space is sparsely covered, PCEvo can effectively enhance the generalization of backbone models.

On MoleculeNet, since the benchmarks are intrinsically small-scale, we evaluate PCEvo under both the standard split (80% training; Table 3) and explicit low-resource splits with training ratios R ∈ 0.1, 0.2, 0.8 (Table 4). The results show that integrating PCEvo into SchNet, DimeNet, and ViSNet leads to consistent gains across tasks and training ratios. The improvement on FreeSolv is especially evident; for example, with SchNet at R=0.2, RMSE decreases from 4.042 to

2.383. Meanwhile, under the standard split, PCEvo also outperforms several widely used pretraining baselines: for instance, PCEvo SchNet achieves an average RMSE of 1.260, improving over MolCLR (1.401) by 0.141.

Overall, PCEvo exhibits stable and transferable improvements across backbones and few-shot settings, and compares favorably to a range of pretraining methods. These results suggest that moving beyond static endpoint regression by decomposing supervision into edit-wise steps and enforcing path-consistency can yield more robust molecular representations in few-shot settings.

4.2 Ablation Study

To isolate the effect of the path-consistency constraint Lcons, we conduct a controlled ablation study. We use SchNet as the backbone, train under the 100-shot, and evaluate on a heldout validation set of 100 disjoint molecules, with ϵHOMO as the target property. Across all runs, we keep the static endpoint supervision Lstatic enabled. Therefore, within each (N, Pmax) configuration, the only difference is whether Lcons is enabled.

Main finding: Lcons is effective when multiple paths are available. Table 5 shows that relative to the single path setting (Pmax = 1), multi-path training (Pmax = 5) consistently improves performance when the path-consistency constraint Lcons is enabled. This trend holds across all neighbor budgets N ∈ 1, 3, 5, indicating that the benefit of Lcons is robust to the choice of N. In contrast, simply increasing the number of paths without Lcons does not lead to reliable gains and can even degrade performance in few-shot settings.

On the role of naive duplication. In our training construction, endpoint molecules are reused across different molecular pairs and paths. Such reuse mainly increases the number of training instances, but does not necessarily introduce additional chemical diversity at the endpoints; therefore, the observed gains cannot be explained by sample count alone. Table 5 provides a controlled comparison supporting this point: once the path-consistency constraint is removed, increasing the path budget (Pmax > 1) or the neighbor budget (N > 1) no longer yields consistent improvements; instead, performance becomes unstable under some configurations.

4.3 Hyperparameter Sensitivity Analysis

To analyze the impact of the maximum number of virtual paths Pmax and the number of source neighbors N on PCEvo, we conduct experiments on the QM9 dataset under the 100shot setting using SchNet as the backbone, with results shown in Fig. 4. Regarding the number of virtual paths, we observe that a single path (Pmax = 1) yields negligible gains over the baseline, indicating that simple augmentation is insufficient. However, increasing Pmax significantly boosts performance by enforcing consistency across diverse topological orderings. This confirms that path-consistent via virtual evolutionary is central to the PCEvo's success. Regarding the context size N, a clear trade-off exists between diversity and relevance. Performance peaks at N = 5; utilizing larger neighborhoods (N ≥ 15) degrades accuracy, as distant neighbors introduce chemically irrelevant scaffolds that generate noisy and overly complex edit paths.

5 Conclusion

In this work, we propose PCEvo, a path-consistent molecular representation method that augments conventional static regression supervision with a virtual evolutionary edit-path consistency constraint. This design enables the model to learn how local structural edits accumulate into global property variations, thereby yielding more robust structure-property mappings. Experiments on QM9 and MoleculeNet under few-shot settings show that PCEvo consistently reduces error, improves stability, and achieves SOTA performance against strong baselines.

Acknowledgments

This work was supported in part by the National Key Research and Development Program of China (2023YFC2705705), the Natural Science Foundation of China (No. 62476203), Key Project of Traditional Chinese Medicine Joint Fund of Hubei Provincial Natural Science Foundation (No.2025AFD47), Hubei Province Science and Technology Innovation Plan Project (No.2025BCB035), the Guangdong Provincial Nat¬ural Science Foundation General Project (No. 2025A1515012155), the Shenzhen Natural Science Foundation Project (No. JCYJ20250604122534006, JCYJ20230807090211021). The Hubei Jiangxia laboratory development project (No. 2025JXKF034),the New Cornerstone Science Foundation through the XPLORER PRIZE,The Innovative Research Group Project of Hubei Province under Grants 2024AFA017.