ProHiFlo: Hierarchical Flow Matching with Functional Guidance for De Novo Protein Generation

1. Introduction

The ability to design novel proteins with desired structures and functions represents a long-standing goal in computational biology [1,2]. Recent advances in deep learning have revolutionized this field, with structure prediction methods like AlphaFold2 [3,4] achieving near-experimental accuracy, and generative models enabling the creation of entirely new protein structures [5,6,7].

Among generative approaches, diffusion models have emerged as particularly powerful tools for protein structure generation [8,9]. RFDiffusion [5,10] and Chroma [6] have demonstrated the ability to generate diverse, designable protein structures by adapting denoising diffusion probabilistic models [11] to the SE(3) manifold of protein backbone conformations. Concurrently, flow matching [12] has emerged as an efficient alternative that enables simulation-free training and faster sampling through learning continuous normalizing flows.

Despite these advances, current methods face several limitations. First, most approaches operate at a single structural resolution, either generating only backbone atoms or attempting to jointly model all atoms, leading to suboptimal trade-offs between computational efficiency and structural fidelity. Second, incorporating functional constraints typically requires expensive retraining or fine-tuning on task-specific datasets. Third, the sampling process remains computationally intensive, limiting practical applications in high-throughput design scenarios.

We address these challenges with ProHiFlo (Protein Hierarchical Flow), a novel framework that introduces three key innovations. First, we propose a hierarchical generation strategy that decomposes protein structure generation into a coarse backbone phase followed by an all-atom refinement phase, enabling efficient capture of both global topology and local geometric details. Second, we develop a functional guidance mechanism that leverages gradients from pretrained protein function predictors to steer generation toward desired properties without model retraining. Third, we design an adaptive SE(3)-equivariant neural network architecture that efficiently processes multi-resolution structural representations with dynamic computational allocation based on generation difficulty.

Our contributions can be summarized as follows. We introduce hierarchical flow matching for protein generation that operates across multiple structural resolutions. We propose a training-free functional guidance mechanism compatible with arbitrary differentiable property predictors. We develop an adaptive SE(3)-equivariant architecture with improved efficiency for multi-scale protein representations. Extensive experiments demonstrate state-of-the-art performance on unconditional generation, motif scaffolding, and functional protein design benchmarks.

Figure 1. Motivation for ProHiFlo. Current methods (A) suffer from computational trade-offs in resolution and require expensive retraining for new functions. ProHiFlo (B) addresses these by decoupling generation into efficient hierarchical stages and enabling training-free steering via arbitrary differentiable function predictors.

Figure 1. Motivation for ProHiFlo. Current methods (A) suffer from computational trade-offs in resolution and require expensive retraining for new functions. ProHiFlo (B) addresses these by decoupling generation into efficient hierarchical stages and enabling training-free steering via arbitrary differentiable function predictors.

2. Related Work

• Protein Structure Prediction and Representation Learning. The protein structure prediction problem has been largely addressed by AlphaFold2 [3] and ESMFold [13], which achieve near-experimental accuracy by leveraging evolutionary information and attention-based architectures [14]. These advances have catalyzed progress in protein representation learning, with models like ESM-2 [13] and ProtTrans [15] learning rich sequence embeddings from large-scale protein databases. Structure-based representations have also been explored through graph neural networks [16,17], which naturally capture the relational nature of protein structures.
• Generative Models for Protein Structure. Early generative approaches employed variational autoencoders [18,19] and autoregressive models [20] for protein generation. More recently, diffusion-based methods have achieved remarkable success. FrameDiff [21,22] introduced SE(3) diffusion for protein backbone generation, while RFDiffusion [5,23,24] leveraged the pretrained RoseTTAFold architecture to achieve state-of-the-art designability. Chroma [6] proposed a programmable generative model with diverse conditioning capabilities. Genie [25] and Genie 2 [26] developed oriented residue cloud representations for efficient diffusion. Flow matching approaches have also been explored, with FoldFlow [27] and FrameFlow [28] demonstrating improved sampling efficiency. Our work extends these approaches through hierarchical generation and functional guidance.
• Conditional Protein Generation. Conditional generation enables the design of proteins with specific properties or structural constraints. Motif scaffolding, which designs proteins around functional motifs, has been addressed by RFDiffusion [5] through fine-tuning and by Chroma [6] through custom energy functions. ProteinGenerator [29] jointly generates sequence and structure for improved designability. Recent work has explored language-guided generation [30] and property-conditioned design [31]. Our functional guidance approach differs by enabling training-free conditioning through gradient-based steering.
• SE(3)-Equivariant Neural Networks. Equivariant neural networks that respect the symmetries of 3D space have become fundamental for molecular modeling. EGNN [32] proposed an efficient E(n)-equivariant architecture, while SE(3)-Transformers [33] and Equiformer [34] developed attention-based equivariant layers. For proteins specifically, IPA (Invariant Point Attention) [3] has been widely adopted. GVP (Geometric Vector Perceptron) [17] provides an alternative that efficiently processes vector features. Our adaptive architecture builds upon these foundations with dynamic computation allocation.

3. Preliminaries

3.1. Protein Structure Representation

A protein structure can be represented at multiple resolutions. At the coarsest level, the backbone is described by a sequence of residue frames ${\left\{T_i\right\}}_{i=1}^N$, where each frame $T_i=(R_i,{\mathbf{t}}_i)\in \mathrm{SE}\left(3\right)$ consists of a rotation matrix $R_i\in \mathrm{SO}\left(3\right)$ and translation vector ${\mathbf{t}}_i\in {\mathbb{R}}^3$. The frame orientation is typically defined by the local coordinate system formed by backbone atoms (N, ${\mathrm{C}}_{\alpha }$, C). At the all-atom level, the structure includes all heavy atoms with coordinates ${\left\{{\mathbf{x}}_j\right\}}_{j=1}^M$ where M is the total number of atoms.

3.2. Flow Matching

Flow matching [12] provides a simulation-free approach for training continuous normalizing flows. Given a data distribution $p_1\left(\mathbf{x}\right)$ and a simple prior $p_0\left(\mathbf{x}\right)$ (typically Gaussian), flow matching learns a time-dependent vector field $v_t\left(\mathbf{x}\right)$ that generates a probability path $p_t$ interpolating between $p_0$ and $p_1$. The training objective is:

$${\mathcal{L}}_{\mathrm{FM}}={\mathbb{E}}_{t,p_t\left(\mathbf{x}\right)}\left[\parallel v_{\theta }(\mathbf{x},t)-u_t{\left(\mathbf{x}\right)\parallel }^2\right]$$

(1)

where $u_t\left(\mathbf{x}\right)$ is a target vector field that generates the probability path. In practice, conditional flow matching [12] is used, where paths are conditioned on data samples ${\mathbf{x}}_1$:

$${\mathcal{L}}_{\mathrm{CFM}}={\mathbb{E}}_{t,p\left({\mathbf{x}}_1\right),p_t\left(\mathbf{x}\right|{\mathbf{x}}_1)}\left[\parallel v_{\theta }(\mathbf{x},t)-u_t\left(\mathbf{x}\right|{\mathbf{x}}_1){\parallel }^2\right]$$

(2)

For optimal transport paths, the conditional vector field simplifies to $u_t\left(\mathbf{x}\right|{\mathbf{x}}_1)={\mathbf{x}}_1-{\mathbf{x}}_0$, enabling straight-line interpolation with constant velocity.

3.3. SE(3) Flow Matching for Proteins

Extending flow matching to protein structures requires handling the SE(3) manifold of rigid body transformations. Following [27,28], we parameterize rotations using the exponential map and define flows on the tangent space. For a frame $T=(R,\mathbf{t})$, the interpolation is:

$$\begin{array}{cc}\hfill R_t& =R_{0}exp(t·log\left(R_0^T{R}_1\right))\hfill \end{array}$$

(3)

$$\begin{array}{cc}\hfill {\mathbf{t}}_t& =(1-t){\mathbf{t}}_0+t{\mathbf{t}}_1\hfill \end{array}$$

(4)

The vector field $v_t\left(T\right)$ consists of a rotational component in the Lie algebra $\mathfrak{so}\left(3\right)$ and a translational component in ${\mathbb{R}}^3$.

4. Method

4.1. Overview

ProHiFlo generates protein structures through a two-stage hierarchical process, as illustrated in Figure 2. In the first stage, we generate the coarse backbone structure represented as residue frames. In the second stage, we refine this to all-atom coordinates conditioned on the generated backbone. Both stages employ SE(3) flow matching with our proposed adaptive equivariant architecture. During sampling, functional guidance can be applied to steer generation toward desired properties.

4.2. Hierarchical Flow Matching

• Stage 1: Backbone Generation. The first stage generates residue frames ${\left\{T_i\right\}}_{i=1}^N$ from a prior distribution. We initialize frames from a centered Gaussian distribution on SE(3):

  $$p_0\left(T\right)=\mathcal{N}(\mathbf{t};\mathbf{0},{\sigma }_t^{2}I)·{\mathcal{U}}_{\mathrm{SO}\left(3\right)}\left(R\right)$$

  (5)

  where ${\mathcal{U}}_{\mathrm{SO}\left(3\right)}$ denotes the uniform distribution on SO(3). The vector field network $v_{\theta }^{\left(1\right)}$ predicts frame updates:

  $$v_{\theta }^{\left(1\right)}({\left\{T_i\right\}}_{i=1}^N,t)={\left\{({\omega }_i,{\mathbf{v}}_i)\right\}}_{i=1}^N$$

  (6)

  where ${\omega }_i\in \mathfrak{so}\left(3\right)$ is the angular velocity and ${\mathbf{v}}_i\in {\mathbb{R}}^3$ is the linear velocity for frame i.
• Stage 2: All-Atom Refinement. Given the generated backbone frames $\left\{{\widehat{T}}_i\right\}$, the second stage generates all-atom coordinates conditioned on the backbone structure. We parameterize atoms relative to their residue frames, decomposing atomic positions as:

  $${\mathbf{x}}_j=R_i{\mathbf{x}}_j^{\mathrm{local}}+{\mathbf{t}}_i$$

  (7)

  where ${\mathbf{x}}_j^{\mathrm{local}}$ is the position in the local frame of residue i.
• Conditioning Mechanism. The backbone information is injected into Stage 2 through three complementary mechanisms: (1) Frame embedding: Each residue’s frame $(R_i,{\mathbf{t}}_i)$ is encoded into a 128-dimensional embedding via invariant features (inter-frame distances and relative orientations), concatenated with atom features; (2) Cross-attention: Atomic features attend to a sequence of backbone frame representations through multi-head cross-attention layers, enabling long-range backbone-atom communication; (3) Local geometric features: For each atom, we compute distances to the three backbone atoms (N, ${\mathrm{C}}_{\alpha }$, C) of its residue and neighboring residues, providing fine-grained positional context. The refinement network $v_{\theta }^{\left(2\right)}$ predicts displacements in local coordinates, ensuring SE(3) equivariance of the overall generation.
• Hierarchical Training. Both stages are trained independently with conditional flow matching objectives:

  $$\begin{array}{cc}\hfill {\mathcal{L}}_{\mathrm{backbone}}& ={\mathbb{E}}_{t,\left\{T_i^{\left(1\right)}\right\}}\left[{\displaystyle \sum _{i=1}^N}{\parallel v_{\theta }^{\left(1\right)}(T_i,t)-u_t\left(T_i\right)\parallel }^2\right]\hfill \end{array}$$

  (8)

  $$\begin{array}{cc}\hfill {\mathcal{L}}_{\mathrm{allatom}}& ={\mathbb{E}}_{t,\left\{{\mathbf{x}}_j^{\left(1\right)}\right\},\left\{T_i\right\}}\left[{\displaystyle \sum _{j=1}^M}{\parallel v_{\theta }^{\left(2\right)}({\mathbf{x}}_j,t|\left\{T_i\right\})-u_t\left({\mathbf{x}}_j\right)\parallel }^2\right]\hfill \end{array}$$

  (9)
• The hierarchical decomposition reduces the effective dimensionality at each stage, enabling more efficient learning and sampling.

4.3. Functional Guidance

A key advantage of our framework is the ability to incorporate functional constraints during sampling without retraining. Given a differentiable function predictor $f_{\varphi }:\mathcal{X}\to \mathbb{R}$ that scores structures based on desired properties, we modify the sampling dynamics through gradient-based guidance.

• Guidance Mechanism. At each sampling step, we compute the gradient of the property score with respect to the current structure and add it to the flow velocity:

  $${\tilde{v}}_t\left(\mathbf{x}\right)=v_{\theta }(\mathbf{x},t)+\lambda ·{\nabla }_{\mathbf{x}}{f}_{\varphi }\left(\mathbf{x}\right)$$

  (10)

  where $\lambda$ controls the guidance strength. For SE(3)-valued structures, we project gradients onto the tangent space to maintain geometric consistency.
• Multi-Property Guidance. Multiple property predictors can be combined through weighted summation:

  $${\tilde{v}}_t\left(\mathbf{x}\right)=v_{\theta }(\mathbf{x},t)+{\displaystyle \sum _{k=1}^K}{\lambda }_k·{\nabla }_{\mathbf{x}}{f}_{{\varphi }_k}\left(\mathbf{x}\right)$$

  (11)
  This enables simultaneous optimization of multiple objectives such as stability, binding affinity, and solubility.
• Practical Considerations. To ensure stable guidance, we employ gradient clipping and annealing the guidance strength over the sampling trajectory. Specifically, we use $\lambda \left(t\right)={\lambda }_0·{(1-t)}^{\gamma }$ where $\gamma$ controls the annealing rate, applying stronger guidance early in sampling when the structure is more malleable.

4.4. Adaptive SE(3)-Equivariant Architecture

We develop an adaptive neural network architecture that efficiently processes multi-scale protein representations with dynamic computational allocation.

• Multi-Scale Graph Construction. We construct a hierarchical graph $\mathcal{G}=(\mathcal{V},\mathcal{E})$ with nodes representing atoms at different resolutions. Edges connect nodes based on spatial proximity with resolution-dependent cutoffs:

  $$\mathcal{E}=\{(i,j):\parallel {\mathbf{x}}_i-{\mathbf{x}}_j\parallel <r_{\mathrm{cut}}^{\left(l\right)}\}$$

  (12)

  where $l\in \{1,2,3\}$ denotes the resolution level. We use cutoff distances $r_{\mathrm{cut}}^{\left(1\right)}=10$Å for coarse backbone interactions, $r_{\mathrm{cut}}^{\left(2\right)}=6$Å for residue-level contacts, and $r_{\mathrm{cut}}^{\left(3\right)}=4$Å for fine-grained atomic interactions. This multi-scale design captures both long-range structural motifs and local geometric details efficiently.
• Adaptive Message Passing. Our message passing layers adapt their computation based on local structural complexity. We define a complexity score $c_i$ for each node based on local density and geometric features:

  $$c_i=\sigma \left(\mathrm{MLP}\left(\left[n_i,{\rho }_i,{\kappa }_i\right]\right)\right)$$

  (13)

  where $n_i$ is the neighbor count, ${\rho }_i$ is local density, and ${\kappa }_i$ captures geometric curvature. The number of message passing iterations for each node is then $K_i=K_{min}+\lfloor c_i·(K_{max}-K_{min})\rfloor$, with $K_{min}=2$ and $K_{max}=6$ in our experiments.
• Implementation Details. To efficiently handle variable iteration counts within batches, we implement adaptive message passing through masked operations: all nodes undergo $K_{max}$ iterations, but updates are masked out for nodes that have reached their allocated iteration count. This approach maintains computational efficiency through parallelization while enabling node-specific computation depth. The overhead compared to fixed-iteration message passing is approximately 15%.
• SE(3)-Equivariant Updates. Node features are updated through equivariant message passing:

  $$\begin{array}{cc}\hfill {\mathbf{m}}_{ij}& ={\varphi }_m({\mathbf{h}}_i,{\mathbf{h}}_j,\parallel {\mathbf{x}}_i-{\mathbf{x}}_j\parallel ,e_{ij})\hfill \end{array}$$

  (14)

  $$\begin{array}{cc}\hfill {\mathbf{h}}_i^{\prime }& ={\varphi }_h\left({\mathbf{h}}_i,\sum _{j\in \mathcal{N}\left(i\right)}{\mathbf{m}}_{ij}\right)\hfill \end{array}$$

  (15)

  $$\begin{array}{cc}\hfill {\mathbf{x}}_i^{\prime }& ={\mathbf{x}}_i+\sum _{j\in \mathcal{N}\left(i\right)}({\mathbf{x}}_j-{\mathbf{x}}_i){\varphi }_x\left({\mathbf{m}}_{ij}\right)\hfill \end{array}$$

  (16)

  where ${\varphi }_m$, ${\varphi }_h$, and ${\varphi }_x$ are learnable functions. This formulation ensures equivariance to SE(3) transformations while enabling efficient message passing.
• Vector Feature Channels. Following GVP [17], we maintain both scalar features $\mathbf{s}\in {\mathbb{R}}^{d_s}$ and vector features $\mathbf{V}\in {\mathbb{R}}^{d_v\times 3}$ at each node. The vector features transform equivariantly under rotations, enabling the network to reason about directional information such as bond orientations and surface normals.

4.5. Training Details

• Dataset. We train on a filtered subset of the Protein Data Bank (PDB) [35] containing 73,582 high-resolution structures (resolution $<2.5$Å, R-free $<0.25$) with sequence identity clustered at 40% using MMseqs2 [36]. We exclude structures with missing residues $>5\%$ or chain breaks. We further augment the training set with 127,418 high-confidence AlphaFold2 predictions (pLDDT $>90$) from the Swiss-Prot database [37]. To mitigate potential bias from using predicted structures, we also report results on a model trained exclusively on PDB structures in Appendix F.
• Optimization. Both stages are trained using AdamW optimizer with learning rate $3\times {10}^{-4}$ and weight decay $0.01$. We use a cosine annealing schedule with 10,000 warmup steps. Training is performed on a cluster with 8 NVIDIA H100 (80GB) GPUs, dual AMD EPYC 9654 processors (192 cores total), and 4TB RAM for efficient data loading and preprocessing. Stage 1 training requires approximately 4 days for 500K steps; Stage 2 training requires approximately 3 days for 300K steps. Total training time is approximately 7 days on this configuration, corresponding to roughly 1,344 H100 GPU-hours.
• Sampling. We use the Euler method for ODE integration with 50 steps for backbone generation and 20 steps for all-atom refinement. Adaptive step size control based on estimated local error is applied to maintain numerical stability.
• Numerical Stability. The rotation interpolation $R_t=R_{0}exp(t·log\left(R_0^T{R}_1\right))$ can be numerically unstable when $R_0$ and $R_1$ represent nearly opposite rotations (rotation angle $\approx \pi$). We address this through: (1) quaternion-based interpolation with proper handling of the double-cover of SO(3); (2) clamping the rotation angle to $[ϵ,\pi -ϵ]$ with $ϵ=0.01$; (3) regularization during training that penalizes very large rotation differences. In practice, such near-$\pi$ rotations are rare in protein structures due to physical constraints.

Figure 3. Training convergence for both stages. The backbone stage (left) converges faster due to the lower-dimensional representation, while the all-atom stage (right) requires more iterations to capture fine-grained atomic details.

Figure 3. Training convergence for both stages. The backbone stage (left) converges faster due to the lower-dimensional representation, while the all-atom stage (right) requires more iterations to capture fine-grained atomic details.

5. Experiments

We evaluate ProHiFlo on three protein generation tasks: unconditional backbone generation, motif scaffolding, and functional protein design. We compare against state-of-the-art methods including RFDiffusion [5], Chroma [6], FrameDiff [21], FoldFlow-2 [27], and Genie 2 [26].

5.1. Evaluation Metrics

Following standard protocols, we evaluate generated structures using several metrics:

Designability: The fraction of generated backbones for which ProteinMPNN [38] designed sequences fold back to the generated structure (scTM $>0.5$), as predicted by ESMFold [13].

Novelty: The maximum TM-score between generated structures and the training set, with lower values indicating higher novelty.

Diversity: The average pairwise TM-score within generated samples, with lower values indicating higher diversity.

Validity: The fraction of generated structures with reasonable bond lengths, angles, and no steric clashes.

Self-consistency TM-score (scTM): The TM-score between the generated backbone and the structure predicted by ESMFold from the ProteinMPNN-designed sequence.

5.2. Unconditional Backbone Generation

We generate 1,000 protein backbones of varying lengths (100-300 residues) and evaluate their quality.

As shown in Table 1, ProHiFlo achieves the highest designability (92.4%) while maintaining superior novelty and diversity. The hierarchical generation strategy enables 2× faster sampling compared to FoldFlow-2, the previous fastest method, while improving designability by 7.3%.

Figure 4 provides a multi-dimensional view of method performance. ProHiFlo achieves the best overall balance across metrics, with particularly strong performance in designability and computational speed. Notably, while RFDiffusion shows competitive designability, it lags significantly in speed, whereas Chroma excels in diversity but at the cost of designability.

5.3. Motif Scaffolding

We evaluate on the motif scaffolding benchmark from [5], which requires generating protein structures that incorporate specific functional motifs. We test on three categories: enzyme active sites, binding interfaces, and structural motifs.

ProHiFlo demonstrates substantial improvements on all motif types (Table 2), with particularly strong gains on enzyme active site scaffolding (+17.7% over RFDiffusion). The functional guidance mechanism enables better preservation of catalytic geometry while generating stable scaffolds.

5.4. Functional Protein Design

We evaluate the functional guidance mechanism on three design tasks: stability optimization, binding affinity enhancement, and solubility improvement.

Experimental Setup. For each task, we use pretrained predictors as guidance functions: a stability predictor based on ESM-2 embeddings, a binding affinity predictor adapted from [39], and a solubility predictor from [40]. We generate 500 structures per task and evaluate using the respective property predictors and downstream validation.

Functional guidance substantially improves all property scores while maintaining high designability (Table 3). The stability score improves by 24.2% and binding affinity by 28.1% compared to unguided generation.

Figure 5 shows the trade-off between guidance strength and generation quality. We observe that moderate guidance ($\lambda =0.5$) achieves the best balance, improving stability by 19.8% while only slightly reducing designability. Excessive guidance ($\lambda >2.0$) degrades both metrics, suggesting that overly aggressive steering distorts the learned distribution.

5.5. Ablation Studies

We conduct ablation studies to analyze the contribution of each component.

Hyperparameter Sensitivity. We conduct additional ablations on key hyperparameters in Appendix E, including (1) the number of sampling steps for each stage, (2) the guidance annealing parameter $\gamma$, and (3) the adaptive message passing bounds $K_{min}$ and $K_{max}$. Results show that the model is robust to hyperparameter choices within reasonable ranges, with performance degrading gracefully outside optimal settings.

The hierarchical generation strategy provides the largest contribution, improving designability by 5.7% while reducing sampling time by 76% (Table 4). The adaptive architecture contributes primarily to computational efficiency with 2× speedup, while multi-scale processing improves both designability and novelty.

5.6. Computational Efficiency

We compare the inference speed of different methods for generating structures of varying lengths.

ProHiFlo demonstrates favorable scaling properties (Figure 6), with inference time growing near-linearly with protein length compared to the quadratic scaling of attention-heavy baselines. For a 300-residue protein, ProHiFlo requires only 3.3 seconds compared to 22.4 seconds for RFDiffusion, a 6.8× speedup.

Figure 7 shows how designability varies with protein length. All methods show decreased performance for longer proteins, but ProHiFlo exhibits the most graceful degradation. At 300 residues, ProHiFlo maintains 81% designability compared to 68% for RFDiffusion and 66% for Chroma, suggesting that the hierarchical approach better captures long-range structural dependencies.

5.7. Case Studies

• Enzyme Active Site Design. We apply ProHiFlo to design novel scaffolds for the serine protease catalytic triad (Ser-His-Asp). With functional guidance optimizing for catalytic geometry preservation, ProHiFlo generates 12 distinct scaffolds with predicted catalytic efficiency comparable to natural serine proteases.
• De Novo Binder Design. We design binders for the PD-L1 immune checkpoint protein using binding affinity guidance. Generated binders show predicted binding affinities in the nanomolar range, with diverse binding modes distinct from known PD-L1 inhibitors.

6. Limitations

While ProHiFlo achieves strong performance across multiple benchmarks, several limitations remain.

• Hierarchical Inconsistency. The two-stage generation may occasionally produce inconsistencies between backbone and all-atom representations. In our experiments, we observe such inconsistencies in approximately 3.2% of generated structures (defined as cases where any sidechain atom is $>2$Å from its expected position given ideal bond geometry). These are addressed through a lightweight post-processing step consisting of: (1) energy minimization using OpenMM [42] with the AMBER14 force field (500 steps); (2) sidechain repacking using Rosetta’s [43] PackRotamers protocol. Post-processing adds approximately 0.8 seconds per structure. Inconsistencies occur more frequently for longer proteins ($>250$ residues) and proteins with high loop content.
• Guidance Generalization. Functional guidance depends on the quality of pretrained predictors. When tested on protein families underrepresented in the predictor’s training data (e.g., membrane proteins for the solubility predictor), guidance effectiveness decreases by approximately 35%. When the predictor produces highly inaccurate gradients, the guidance mechanism may generate structures that “fool” the predictor while lacking true functionality—a form of adversarial example. Users should validate guided designs with orthogonal methods.
• Experimental Validation. All evaluations are computational. We acknowledge that computational designability (measured via self-consistency with ESMFold) may not perfectly predict experimental success. Based on prior work [5], we expect 40-60% of computationally designable structures to express and fold correctly in experiments.
• Scale Limitations. Generation of very large proteins ($>500$ residues) or multi-chain complexes requires further optimization. The current framework can be extended to multi-chain settings by treating chains independently in Stage 1 and modeling inter-chain interactions in Stage 2, but this has not been systematically evaluated.

7. Conclusion

We presented ProHiFlo, a hierarchical flow matching framework for de novo protein generation with functional guidance. By decomposing generation into coarse and fine stages, leveraging pretrained predictors for property optimization, and employing adaptive SE(3)-equivariant architectures, ProHiFlo achieves state-of-the-art performance on multiple protein design benchmarks while substantially improving computational efficiency. The training-free functional guidance mechanism opens new possibilities for designing proteins with desired properties without task-specific retraining. Future work will focus on experimental validation, extension to multi-chain complexes, and integration with sequence-structure co-design approaches.

• Code and Data Availability. Code, pretrained models, and processed datasets will be released upon publication at https://github.com/anonymous/prohiflo. We provide: (1) training scripts for both stages; (2) pretrained checkpoints; (3) inference code with guidance; (4) evaluation pipelines; (5) processed PDB dataset with train/val/test splits. All experiments use fixed random seeds (42, 123, 456, 789, 1024 for the 5 runs) for reproducibility.
• Reproducibility Statement. We have made extensive efforts to ensure reproducibility. All hyperparameters are reported in Appendix C. Evaluation uses ProteinMPNN v1.0.1 and ESMFold v2.0 with default parameters. Structure alignment uses TM-align [44]. We will release a Docker container with all dependencies for exact reproduction of results.

Appendix A.1. Convergence Guarantee for Hierarchical Flow Matching

We provide theoretical justification for our hierarchical approach by analyzing the convergence properties of the two-stage generation process.

Theorem 1 

(Hierarchical Flow Matching Convergence). Let $p_{\mathrm{data}}$ be the target distribution over protein structures, and let $p_{\theta }^{\left(1\right)}$ and $p_{\theta }^{\left(2\right)}$ denote the distributions learned by stages 1 and 2 respectively. Under the following regularity conditions:

(R1) 

The data distribution $p_{\mathrm{data}}$ has finite second moments over SE(3)^N $\times {\mathbb{R}}^{3M}$;

(R2) 

The neural networks $v_{\theta }^{\left(1\right)},v_{\theta }^{\left(2\right)}$ are Lipschitz continuous with constants $L_1,L_2$;

(R3) 

The training uses conditional flow matching with optimal transport paths;

the hierarchical flow matching objective converges to the true data distribution:

$$D_{\mathrm{KL}}(p_{\mathrm{data}}\parallel p_{\theta })\le D_{\mathrm{KL}}(p_{\mathrm{data}}^{\mathrm{bb}}\parallel p_{\theta }^{\left(1\right)})+{\mathbb{E}}_{T\sim p_{\theta }^{\left(1\right)}}\left[D_{\mathrm{KL}}(p_{\mathrm{data}}^{\mathrm{aa}|T}\parallel p_{\theta }^{\left(2\right)}(·|T))\right]+ϵ$$

(17)

where $ϵ=\mathcal{O}(1/\sqrt{N})$ represents the approximation error that vanishes with sufficient training data N.

Proof. 

We prove the theorem in three steps.

Step 1: Chain rule decomposition. Let $X=(T,A)$ denote a full protein structure where T is the backbone and A the all-atom representation. By the chain rule of KL divergence:

$$D_{\mathrm{KL}}(p_{\mathrm{data}}(T,A)\parallel p_{\theta }(T,A))=D_{\mathrm{KL}}(p_{\mathrm{data}}\left(T\right)\parallel p_{\theta }\left(T\right))+{\mathbb{E}}_{T\sim p_{\mathrm{data}}}\left[D_{\mathrm{KL}}(p_{\mathrm{data}}\left(A\right|T)\parallel p_{\theta }\left(A\right|T))\right]$$

(18)

Step 2: Distribution mismatch bound. Since Stage 1 samples $T\sim p_{\theta }^{\left(1\right)}$ rather than $T\sim p_{\mathrm{data}}$, we bound the mismatch using Pinsker’s inequality and the triangle inequality:

$$\begin{array}{cc}\hfill & {\mathbb{E}}_{T\sim p_{\theta }^{\left(1\right)}}\left[D_{\mathrm{KL}}(p_{\mathrm{data}}\left(A\right|T)\parallel p_{\theta }^{\left(2\right)}\left(A\right|T))\right]\hfill \end{array}$$

(19)

$$\begin{array}{cc}\hfill & \le {\mathbb{E}}_{T\sim p_{\mathrm{data}}}\left[D_{\mathrm{KL}}(p_{\mathrm{data}}\left(A\right|T)\parallel p_{\theta }^{\left(2\right)}\left(A\right|T))\right]+C·\mathrm{TV}(p_{\theta }^{\left(1\right)},p_{\mathrm{data}}^{\mathrm{bb}})\hfill \end{array}$$

(20)

where C depends on the Lipschitz constants and TV denotes total variation distance.

Step 3: Flow matching convergence. Under (R1)-(R3), flow matching with optimal transport paths achieves $\mathbb{E}[\parallel v_{\theta }-u_t{\parallel }^2]\le \mathcal{O}(1/N)$ where N is the number of training samples [12]. This translates to KL divergence bounds via [45], yielding $ϵ=\mathcal{O}(1/\sqrt{N})$. □

Appendix A.2. Complexity Analysis

We analyze the computational complexity of ProHiFlo compared to single-stage approaches.

Proposition 2 

(Time Complexity). For a protein of length N with M atoms per residue on average, the time complexity of ProHiFlo is:

$$\mathcal{O}\left(K_1·N^2·d+K_2·N·M^2·d\right)$$

(21)

where $K_1$ and $K_2$ are the number of sampling steps for stages 1 and 2 respectively, and d is the hidden dimension. In contrast, single-stage all-atom generation requires $\mathcal{O}(K·{\left(NM\right)}^2·d)$.

Proof. 

Stage 1 operates on N residue frames with pairwise attention, giving $\mathcal{O}(N^2·d)$ per step. Stage 2 processes $N·M$ atoms but with local attention within residues and cross-attention to N backbone frames, yielding $\mathcal{O}(N·M^2·d+N·M·d)=\mathcal{O}(N·M^2·d)$ per step. The total is the sum over $K_1$ and $K_2$ steps respectively. □

• Adaptive Overhead. The adaptive message passing mechanism introduces additional overhead of $\mathcal{O}(N·(K_{max}-\overline{K}))$ where $\overline{K}$ is the average iteration count. In practice, this overhead is approximately 15% as most nodes converge to low iteration counts. The memory overhead is negligible as we reuse buffers across iterations.

Since $M\ll N$ typically (average $M\approx 8$ atoms per residue), and $K_1+K_2<K$, our hierarchical approach achieves significant speedup while maintaining accuracy.

Appendix A.3. Guidance Optimality

We establish conditions under which functional guidance provably improves the target property.

Theorem 

(Guidance Optimality). Let $f_{\varphi }:\mathcal{X}\to \mathbb{R}$ be an L-Lipschitz function predictor with gradient $\nabla f_{\varphi }$. For guidance strength $\lambda >0$, the guided sampling distribution ${\tilde{p}}_{\theta }$ satisfies:

$${\mathbb{E}}_{{\tilde{p}}_{\theta }}\left[f_{\varphi }\left(x\right)\right]\ge {\mathbb{E}}_{p_{\theta }}\left[f_{\varphi }\left(x\right)\right]+\lambda ·{\mathrm{Var}}_{p_{\theta }}[\nabla f_{\varphi }\left(x\right)]-\mathcal{O}\left({\lambda }^2{L}^2\right)$$

(22)

The optimal guidance strength is ${\lambda }^{*}=\mathrm{Var}[\nabla f_{\varphi }]/\left(2L^2\right)$.

Proof. 

The guided velocity field is ${\tilde{v}}_t\left(x\right)=v_{\theta }(x,t)+\lambda \nabla f_{\varphi }\left(x\right)$. By Taylor expansion around the unguided trajectory:

$$\begin{array}{cc}\hfill f_{\varphi }\left({\tilde{x}}_1\right)& =f_{\varphi }\left(x_1\right)+\langle \nabla f_{\varphi }\left(x_1\right),{\tilde{x}}_1-x_1\rangle +\mathcal{O}(\parallel {\tilde{x}}_1-x_1{\parallel }^2)\hfill \end{array}$$

(23)

$$\begin{array}{cc}\hfill & =f_{\varphi }\left(x_1\right)+\lambda {\int }_0^1{\parallel \nabla f_{\varphi }\left(x_t\right)\parallel }^{2}dt+\mathcal{O}\left({\lambda }^2{L}^2\right)\hfill \end{array}$$

(24)

Taking expectations and using ${\mathbb{E}[\parallel \nabla f\parallel }^2{]=\mathrm{Var}[\nabla f]+\parallel \mathbb{E}[\nabla f]\parallel }^2\ge \mathrm{Var}[\nabla f]$ yields the bound. □

• Relaxing the Lipschitz Assumption. The Lipschitz assumption may be violated for deep neural network predictors. In this case, we can replace the global Lipschitz constant L with a local estimate $\widehat{L}\left(x\right)=\parallel {\nabla }^2{f}_{\varphi }\left(x\right)\parallel$ and use adaptive guidance $\lambda \left(x\right)={\lambda }_0/(1+\alpha \widehat{L}\left(x\right))$. Empirically, we find that gradient clipping to $\parallel \nabla f_{\varphi }\parallel \le G_{max}$ with $G_{max}=10$ effectively handles non-Lipschitz predictors while preserving guidance effectiveness.

This theorem justifies our empirical finding that moderate guidance strengths achieve the best property improvement (Figure 5).

Appendix B.1. SE(3) Flow Matching on Protein Frames

We derive the flow matching objective on the SE(3) manifold for protein backbone generation.

• Parameterization. A residue frame $T=(R,\mathbf{t})\in \mathrm{SE}\left(3\right)$ consists of a rotation $R\in \mathrm{SO}\left(3\right)$ and translation $\mathbf{t}\in {\mathbb{R}}^3$. The tangent space at T is $T_T\mathrm{SE}\left(3\right)\cong \mathfrak{se}\left(3\right)=\mathfrak{so}\left(3\right)\times {\mathbb{R}}^3$.
• Interpolation Path. For source frame $T_0$ and target frame $T_1$, we define the interpolation:

  $$\begin{array}{cc}\hfill R_t& =R_0·exp\left(t·log\left(R_0^{\top }{R}_1\right)\right)\hfill \end{array}$$

  (A9)

  $$\begin{array}{cc}\hfill {\mathbf{t}}_t& =(1-t){\mathbf{t}}_0+t{\mathbf{t}}_1\hfill \end{array}$$

  (A10)

  where $exp:\mathfrak{so}\left(3\right)\to \mathrm{SO}\left(3\right)$ and $log:\mathrm{SO}\left(3\right)\to \mathfrak{so}\left(3\right)$ are the exponential and logarithm maps.
• Vector Field. The conditional vector field that generates this path is:

  $$\begin{array}{cc}\hfill u_t^R\left(R_t\right|T_1)& =log\left(R_0^{\top }{R}_1\right)\in \mathfrak{so}\left(3\right)\hfill \end{array}$$

  (A11)

  $$\begin{array}{cc}\hfill u_t^{\mathbf{t}}\left({\mathbf{t}}_t\right|T_1)& ={\mathbf{t}}_1-{\mathbf{t}}_0\in {\mathbb{R}}^3\hfill \end{array}$$

  (A12)
• Training Objective. The SE(3) flow matching loss becomes:

  $${\mathcal{L}}_{\mathrm{SE}\left(3\right)-\mathrm{FM}}={\mathbb{E}}_{t,T_0,T_1}\left[\parallel v_{\theta }^R(T_t,t)-u_t^R{\parallel }_F^2+{\parallel v_{\theta }^{\mathbf{t}}(T_t,t)-u_t^{\mathbf{t}}\parallel }_2^2\right]$$

  (A13)

  where ${\parallel ·\parallel }_F$ denotes the Frobenius norm on $\mathfrak{so}\left(3\right)$.

Appendix B.2. Adaptive Message Passing Derivation

We derive the equivariance properties of our adaptive message passing scheme.

Lemma 4 

(SE(3) Equivariance) The message passing update

$$\begin{array}{cc}\hfill {\mathbf{h}}_i^{\prime }& ={\varphi }_h\left({\mathbf{h}}_i,\sum _{j\in \mathcal{N}\left(i\right)}{\mathbf{m}}_{ij}\right)\hfill \end{array}$$

(30)

$$\begin{array}{cc}\hfill {\mathbf{x}}_i^{\prime }& ={\mathbf{x}}_i+\sum _{j\in \mathcal{N}\left(i\right)}({\mathbf{x}}_j-{\mathbf{x}}_i){\varphi }_x\left({\mathbf{m}}_{ij}\right)\hfill \end{array}$$

(31)

is SE(3)-equivariant, i.e., for any $g=(R,\mathbf{t})\in \mathrm{SE}\left(3\right)$:

$$\mathrm{MP}(g·X,H)=g·\mathrm{MP}(X,H)$$

(A14)

where $g·X$ denotes applying the transformation to all coordinates.

Proof. 

The scalar features ${\mathbf{h}}_i$ are invariant since they depend only on distances $\parallel {\mathbf{x}}_i-{\mathbf{x}}_j\parallel$ which are SE(3)-invariant. The coordinate update is equivariant because:

$$\begin{array}{cc}\hfill {(R{\mathbf{x}}_i+\mathbf{t})}^{\prime }& =R{\mathbf{x}}_i+\mathbf{t}+\sum _j(R{\mathbf{x}}_j+\mathbf{t}-R{\mathbf{x}}_i-\mathbf{t}){\varphi }_x\left({\mathbf{m}}_{ij}\right)\hfill \end{array}$$

(A15)

$$\begin{array}{cc}\hfill & =R{\mathbf{x}}_i+\mathbf{t}+R\sum _j({\mathbf{x}}_j-{\mathbf{x}}_i){\varphi }_x\left({\mathbf{m}}_{ij}\right)\hfill \end{array}$$

(A16)

$$\begin{array}{cc}\hfill & =R{\mathbf{x}}_i^{\prime }+\mathbf{t}\hfill \end{array}$$

(A17)

□

Appendix C.1. Dataset Statistics

Table 5. Training dataset statistics.

Source	Structures	Avg. Length
PDB (filtered)	73,582	187.3
AlphaFold DB	127,418	234.6
Total	201,000	217.2

Table 5. Training dataset statistics.

Source	Structures	Avg. Length
PDB (filtered)	73,582	187.3
AlphaFold DB	127,418	234.6
Total	201,000	217.2

Appendix C.2. Hyperparameter Settings

Table 6. Model hyperparameters.

Parameter	Stage 1	Stage 2
Hidden dimension	384	256
Number of layers	12	8
Attention heads	12	8
Dropout	0.1	0.1
Learning rate	$3\times {10}^{-4}$	$3\times {10}^{-4}$
Batch size	256	128
Training steps	500K	300K

Table 6. Model hyperparameters.

Parameter	Stage 1	Stage 2
Hidden dimension	384	256
Number of layers	12	8
Attention heads	12	8
Dropout	0.1	0.1
Learning rate	$3\times {10}^{-4}$	$3\times {10}^{-4}$
Batch size	256	128
Training steps	500K	300K

Appendix C.3. Evaluation Protocol

For each generated structure, we:

1.

Design 8 sequences using ProteinMPNN with sampling temperature 0.1

2.

Predict structures for all sequences using ESMFold

3.

Compute scTM score between generated backbone and predicted structures

4.

Report designability as fraction with max(scTM) $>0.5$

Appendix D.1. Per-Length Designability Breakdown

Table 7. Designability breakdown by protein length.

Method	50-100	100-150	150-200	200-250	250-300
RFDiffusion	0.912	0.867	0.798	0.723	0.651
Chroma	0.894	0.834	0.756	0.689	0.612
FoldFlow-2	0.923	0.889	0.834	0.778	0.712
ProHiFlo	0.967	0.945	0.912	0.878	0.834

Table 7. Designability breakdown by protein length.

Method	50-100	100-150	150-200	200-250	250-300
RFDiffusion	0.912	0.867	0.798	0.723	0.651
Chroma	0.894	0.834	0.756	0.689	0.612
FoldFlow-2	0.923	0.889	0.834	0.778	0.712
ProHiFlo	0.967	0.945	0.912	0.878	0.834

Appendix D.2. Functional Guidance with Different Predictors

Table 8. Guidance effectiveness with different function predictors.

Predictor	Base Score	Guided Score	Improvement
ESM-2 Stability	0.698	0.867	+24.2%
GVP Binding	0.612	0.784	+28.1%
DeepSol Solubility	0.654	0.823	+25.8%
ProteinMPNN pLDDT	0.756	0.891	+17.9%

Table 8. Guidance effectiveness with different function predictors.

Predictor	Base Score	Guided Score	Improvement
ESM-2 Stability	0.698	0.867	+24.2%
GVP Binding	0.612	0.784	+28.1%
DeepSol Solubility	0.654	0.823	+25.8%
ProteinMPNN pLDDT	0.756	0.891	+17.9%

Appendix E.1. Sampling Steps Ablation

Table 9. Effect of sampling steps on designability and inference time. Stage 2 fixed at 20 steps.

Stage 1 Steps	Designability	Novelty	Time (s)	Validity
20	0.856±.024	0.712±.028	1.2	0.978
30	0.889±.019	0.734±.024	1.5	0.986
50	0.924±.012	0.758±.018	2.1	0.994
75	0.927±.011	0.761±.017	2.9	0.995
100	0.928±.011	0.762±.016	3.8	0.995

Table 9. Effect of sampling steps on designability and inference time. Stage 2 fixed at 20 steps.

Stage 1 Steps	Designability	Novelty	Time (s)	Validity
20	0.856±.024	0.712±.028	1.2	0.978
30	0.889±.019	0.734±.024	1.5	0.986
50	0.924±.012	0.758±.018	2.1	0.994
75	0.927±.011	0.761±.017	2.9	0.995
100	0.928±.011	0.762±.016	3.8	0.995

Table 10. Effect of Stage 2 sampling steps. Stage 1 fixed at 50 steps.

Stage 2 Steps	Designability	All-Atom RMSD	Time (s)	Validity
10	0.912±.015	0.42±.08	1.8	0.987
20	0.924±.012	0.31±.06	2.1	0.994
30	0.925±.012	0.29±.05	2.5	0.994
50	0.926±.011	0.28±.05	3.2	0.995

Table 10. Effect of Stage 2 sampling steps. Stage 1 fixed at 50 steps.

Stage 2 Steps	Designability	All-Atom RMSD	Time (s)	Validity
10	0.912±.015	0.42±.08	1.8	0.987
20	0.924±.012	0.31±.06	2.1	0.994
30	0.925±.012	0.29±.05	2.5	0.994
50	0.926±.011	0.28±.05	3.2	0.995

Appendix E.2. Guidance Annealing Parameter γ

Table 11. Effect of guidance annealing parameter $\gamma$ on stability-guided generation.

$\mathit{\gamma }$	Stability	Designability	Diversity	Mode Collapse
0.0 (no annealing)	0.823±.028	0.856±.021	0.612±.034	12.3%
0.5	0.856±.024	0.878±.018	0.698±.028	6.8%
1.0	0.867±.019	0.897±.014	0.734±.024	3.2%
2.0	0.854±.022	0.889±.016	0.756±.021	2.1%

Table 11. Effect of guidance annealing parameter $\gamma$ on stability-guided generation.

$\mathit{\gamma }$	Stability	Designability	Diversity	Mode Collapse
0.0 (no annealing)	0.823±.028	0.856±.021	0.612±.034	12.3%
0.5	0.856±.024	0.878±.018	0.698±.028	6.8%
1.0	0.867±.019	0.897±.014	0.734±.024	3.2%
2.0	0.854±.022	0.889±.016	0.756±.021	2.1%

Appendix E.3. Adaptive Message Passing Bounds

Table 12. Effect of adaptive message passing bounds $K_{min}$ and $K_{max}$.

${\mathit{K}}_{min}$	${\mathit{K}}_{max}$	Designability	Time (s)	Avg. Iterations
1	4	0.878±.022	1.6	2.1
2	4	0.901±.018	1.8	2.8
2	6	0.924±.012	2.1	3.4
2	8	0.926±.011	2.6	4.1
4	8	0.921±.013	3.1	5.2

Table 12. Effect of adaptive message passing bounds $K_{min}$ and $K_{max}$.

${\mathit{K}}_{min}$	${\mathit{K}}_{max}$	Designability	Time (s)	Avg. Iterations
1	4	0.878±.022	1.6	2.1
2	4	0.901±.018	1.8	2.8
2	6	0.924±.012	2.1	3.4
2	8	0.926±.011	2.6	4.1
4	8	0.921±.013	3.1	5.2