Random Feature Representation Boosting

Thank you for taking the time to review our work and for your comments. Below, we address your questions and clarify the novelty and key differences of our work.

Q: What are the main differences between Random Feature Representation Boosting and existing Functional Boosting? Are dense Random Feature Representation Boosting and Functional Gradient Boosting the same?

We appreciate the opportunity to clarify the distinctions. There are crucial differences between RFRBoost and what might be considered 'Existing Functional Gradient Boosting' (FGB), particularly the classical form used in methods like XGBoost.

• Boosting Target (Label Space vs. Feature Space): Classical FGB, as outlined in our Section 2.1, typically operates in label space. It iteratively adds weak learners (like decision trees) to directly improve the model's prediction $F_t(x)$ by fitting the gradient of the loss with respect to the current prediction. In contrast, RFRBoost operates in the framework of Gradient Representation Boosting (GRB), detailed in Section 2.2. GRB boosts in feature space, iteratively building a deep residual-neural-network-like feature representation $\Phi_t(x)$ by adding residual blocks $g_t$ that approximate the functional gradient of the loss with respect to the neural network representation itself. The final prediction is then made by a single linear model $W^T \Phi_T(x)$ trained on the final representation. RFRBoost in the $\Delta_{dense}$ regime is not equivalent to classical gradient boosting, because only a single linear model is trained on the boosted ResNet feature representation, rather than an ensemble of such models as in traditional gradient boosting.

• Random Feature Residual Blocks $g_t$: While the theory of GRB allows for general function classes $g_t$, RFRBoost specifically uses fixed random features, allowing for optimized computations and specialized theoretical results. More specifically, we have $g_t = A_t f_t$, where $f_t = f_t(x, \Phi_{t-1}(x))$ is generated using random weights, and only the linear map $A_t$ is trained using linear methods. This is fundamentally different from previous GRB methods that train bigger neural networks end-to-end as residual blocks using SGD. Our use of fixed random features is what enables efficient computational routines and allows for analytical solutions as per Theorems 3.1 and 3.2, while retaining the theoretical guarantees stemming from GRB (Theorem 3.3), and significantly improving performance (Section 4.1). A specific contribution within our gradient-greedy approach is the incorporation of a unit-norm constraint when learning $g_t$, differing from related literature. This ensures we correctly identify the optimal direction for the residual block in function space, which is crucial for the validity of the functional Taylor approximation of the risk. The special setting of random features is what allowed us to prove that this can be obtained by solving a constrained least squares problem (Theorem 3.2).

Q: Can you provide insights into why XGBoost achieves a lower RMSE in regression tasks?

It is challenging to pinpoint a single definitive reason why XGBoost achieves slightly lower RMSE than RFRBoost in our experiments, as RFRBoost and XGBoost operate under fundamentally different paradigms as explained above. (Note however that while XGBoost achieves a slightly lower RMSE than RFRBoost, it performs worse in terms of average rank.)

XGBoost uses classical gradient boosting (Section 2.1) to construct a large ensemble, comprising hundreds or thousands of decision trees. It learns by fitting trees to the residual errors (gradients) in label/logit space. Its strengths lie in its ability to capture complex interactions via hierarchical, axis-aligned splits and regularization techniques tailored for trees.

RFRBoost, conversely, uses gradient representation boosting (Section 2.2) to build a deep feature representation $\Phi_T(x)$ using random projections, followed by a single linear predictor. Its inductive bias stems from the nature of the random features used (e.g., dense projections with tanh activation in our experiments) and the layer-wise boosting process in feature space.

Potential factors contributing to the observed difference on these specific benchmarks might be due to different inductive biases of the different approaches, ensembling vs using a single predictor, or difference in regularization of discrete trees versus continuous (logistic) regression targets.

---

We hope this clarifies the methodological distinctions and provides context for the performance comparison. We aim to make this clearer in the revised paper.

We sincerely thank you for your thorough and positive assessment of our work and for your valuable comments. Below, we address your questions:

Empirical validation primarily focuses on tabular data: The main question we set out to answer in our paper was whether one could build deep random feature neural networks (RFNNs) that significantly improve the performance of single-layer RFNNs while retaining their computational benefits. Since traditional RFNNs, using dense random projections, are primarily designed for tabular data, our paper focused on this domain. A careful study on tabular data is a necessary first step before looking at structured data like images or sequences, that would necessitate more complex random feature structures (e.g., random convolutional kernels or recurrent structures). Integrating such major architectural changes introduces complexities that make it difficult to address thoroughly within a single theory and within the scope of this single paper. We therefore consider this a valuable avenue for future research.

Q: The paper hypothesizes briefly that this improvement may arise due to the constraint on the functional gradient norm, but it lacks deeper analytical or empirical exploration of why this occurs.

We appreciate the request for further clarification, however we are not quite sure what "improvement" here is referring to. We give two possible interpretations:

1: Regarding the empirical observation in Table 1 where the gradient-greedy approach slightly outperforms the exact-greedy approach for regression tasks, we hypothesise that this could be due to an implicit regularisation effect. The gradient-greedy method, by focusing only on aligning with the functional gradient direction (under a norm constraint), might avoid overfitting compared to the exact-greedy method which directly minimizes the loss at each step. Providing a rigorous theoretical analysis of this phenomenon is challenging and would likely require fundamentally new tools, making it an interesting direction for future theoretical work.

2: Regarding why our gradient-greedy approach succeeds while Suggala et al. (2020) observed it performed worse for medium-sized SGD-trained networks, we refer to our discussion at the end of Section 2.2. The crucial difference lies in how the gradient step is implemented. Our gradient-greedy approach (Algorithm 2 and Theorem 3.2) explicitly solves an optimization problem constrained to find the unit-norm functional direction of the residual block $g$ that best aligns with the negative functional gradient $-\nabla \widehat{R}$ (see Theorem 3.2). This ensures that the learned residual block $g$ points in the correct direction in function space, allowing for a valid functional first-order approximation used in boosting (Eq. 4). In contrast, standard SGD training of a residual block to minimize the functional inner product $\langle g, \nabla_2 \widehat{R} \rangle_{L_2^D({\mu})}$ does not inherently enforce a constraint on the functional norm. We believe incorporating an explicit norm constraint, focusing on finding the optimal direction, is key to the success of the gradient-greedy strategy in our random feature setting.

Dataset Size and Scalability: Please see `Dataset Size and Scalability' in response to Reviewer XEnQ.

Originality and Contribution We respectfully disagree that our contribution is merely an integration of established components. While random features, residual blocks, and boosting are known concepts, successfully constructing deep RFNNs presents unique and significant challenges, since naively stacking random layers has been observed to degrade performance. Our key contribution is providing the first principled framework for building deep RFNNs using any general random feature type by rigorously using ideas from gradient representation boosting (which differs from classical gradient boosting). The use of random features within this framework is not just an implementation detail; it is fundamental to enabling the computational efficiency and theoretical analysis central to RFRBoost. This includes the closed-form solutions of Theorem 3.1, and the theory leading to the quadratically constrained least squares problem of Theorem 3.2.

Reviewer Comment: No supplementary material: We would like to clarify that supplementary material containing the Python code for our model implementations (RFRBoost, baselines) and experimental procedures was provided with our submission. We commit to making this code publicly available in a repository upon acceptance to ensure full reproducibility.

---

We hope these clarifications adequately address your comments and further highlight the contributions of our work.

We appreciate your valuable comments and thoughtful points raised, especially the positive feedback regarding experimental rigour. We address your points below.

Are results MSE only?: We apologize if this was unclear. Regression tasks are reported in Table 1 and Figure 2, while classification tasks (cross-entropy loss) are reported in Table 2 and Figure 3. We will revise the text in Section 4.1 to explicitly state that the regression and classification results are presented separately. Additionally, we will also add scatter/box-and-whiskers plots to better display the variability of the models, an ablation study for SWIM vs iid features, and full dataset-wise results in the Appendix.

RFRBoost supports any convex loss function, and the main body of the paper was developed with this in mind. The exact-greedy analytical solution are a special case for MSE loss, but the gradient-greedy algorithm of Sec 3.3 and Algo. 2 supports general convex loss function. More importantly, in the gradient-greedy case we argued why previous work on SGD-based representation boosting overlooked the crucial detail of constraining the functional norm of the residual block to unit norm. This is essential to obtain a valid functional Taylor approximation of the risk, and we proved in Theorem 3.2 how to efficiently calculate and incorporate this constraint when using random features. This step provably involves solving a constrained least squares problem (where the functional gradient vector depends on the convex loss), which is where the MSE confusion might have stemmed from.

References on CNNs and scattering networks: We thank the referee for pointing out this related field of study. We will incorporate these topics into the related works section in the revised manuscript. We agree that there is no immediate way to include these in experimentation as they concern image data. However, this is an interesting and important direction for future research to build upon the foundational work presented in our paper.

Code Blocks: We are glad you found the embedded Python code snippets helpful for clarity and reproducibility. Following your suggestion, we have added corresponding code examples for the gradient calculation in Appendix C.

Dataset Size and Scalability: The decision to truncate datasets at 5000 observations stemmed primarily from computational resource constraints, balanced against the desire to experiment on a large number (91) of diverse datasets. Note that this threshold is substantially larger than the median dataset size (approx. 2400) within the OpenML Benchmark Suite. The primary computational bottlenecks were the E2E ResNet, and the rigorous nested 5-fold cross-validation procedure with Optuna (requiring 5*100 fits per dataset fold per model). Conducting such extensive evaluation on 91 datasets at their full size was computationally infeasible for us. However, we argue that evaluating performance on datasets of this scale is valuable, as this regime is highly relevant for many practical applications and is a common domain for RFNNs and extreme learning machines. Restricting dataset size is not uncommon in published papers (e.g. Bolager et al. 2023); we aimed for transparency by stating it directly in the main paper, rather than leaving this to the appendix.

To further address scalability we will:

• Add an analysis of the computational complexity of Algorithms 1 and 2, noting its scaling behaviour $O(NTD^2p^2d)$ for training.

• Add experiments in the Appendix where RFRBoost and baselines are trained on larger datasets using increasing fractions of the data (e.g., 10k, 25k, 50k, 100k, 200k, 400k samples, Covertype, YP-MSD). These will be displayed in a plot with number of samples on the X axis, and results/time on Y, to demonstrate scaling behaviour.

Scaling RFRBoost efficiently to datasets with millions of rows would likely require distributed computing paradigms such as MapReduce and implementation optimizations beyond the scope of this initial work, which we identify as an important direction for future research.

Terminology: We have replaced the phrase "law of the data" with "distribution of the data". The original phrasing, while common in some areas of theoretical statistics, might be less familiar to the broader machine learning audience, and we have revised the manuscript accordingly.

Sensitivity of Parameters: We observed that optimal hyperparameters were highly dataset-dependent for all models, including XGBoost, ridge/logistic regression, E2E MLP ResNets, and RFRBoost. This variability motivated our use of a rigorous nested 5-fold CV scheme coupled with Optuna for hyperparameter tuning, ensuring a fair comparison on all datasets. Anecdotally, we observed a general trend where RFRBoost performance tended to improve with increasing depth, whereas the optimal E2E MLP ResNets depth was more varied.

---

We hope these responses adequately address your questions.

We sincerely thank Reviewer fxoT for their careful reading and constructive feedback. Below, we provide detailed responses to each of the points raised.

On Critical Difference Diagrams: To facilitate a meaningful comparison across all baseline models and RFRBoost variants, we use the Wilcoxon signed-rank test with Holm correction, as mentioned in Sec. 4.1 Evaluation Procedure. This statistical test assesses significant pairwise differences between models, and the average relative ranks and statistically significant groupings are presented visually in critical difference diagrams (Figures 2 and 3). If two models are connected by bars, then the difference between the two methods is not statistically significant. This is a well-established methodology for comparing multiple algorithms across multiple datasets [Demšar, 2006; Benavoli et al., 2016] and is standard practice within the broader statistical learning and data mining literature (see also the enthusiasm of Reviewer JjfD regarding the use of Wilcoxon signed-rank test). We aim to explain this in more detail in the revised manuscript.

Presentation of Results Tables: We agree that presenting only mean scores limits the interpretation of the results. To address this, we will incorporated your suggestion of creating scatter plots overlaid with box-and-whiskers plots to visualize the distribution and variability of performance across all 91 datasets for each method. The critical difference diagrams (Figures 2 and 3) complement this by indicating which observed differences in relative rank are statistically significant. Furthermore, we will add a section to the Appendix with complete dataset-wise results for each model. All regression targets have been normalized to have 0 mean and variance 1 (Section 4.1), hence the results should be relatively comparable.

Dataset Size & Scalability: Please see `Dataset Size and Scalability' in the response to Reviewer XEnQ.

SWIM initialization vs iid: We performed experiments using both iid random features and SWIM initialized random features. We initially presented only the SWIM results in the main paper due to slightly superior empirical performance and the 8 page space constraint. In the revised manuscript, we will add an ablation study with pairwise plots directly comparing RFRBoost with iid features versus SWIM features, quantifying the (small) difference in performance. For example, SWIM features only beat iid features 51 out of 91 times.

How were the methods implemented? We will add this missing detail to the revised manuscript. All baseline models were implemented in PyTorch, and RFRBoost in particular was implemented following Algorithms 1 and 2. We used LBFG-S as the minimizer for log likelihood loss, and standard torch linear algebra routines for the (sandwiched) least squares, following our derived formulas in the Appendix. E2E ResNets were trained with Adam as described in Section 4.1. The full code is currently contained in the supplementary material and we will additionally release it in a public repository upon acceptance.

Pseudoinverse: Thanks for this remark. We will point out the alternative expression in terms of the pseudoinverse in the revised version.

Theoretical Regret Bound: Our regret bound in Theorem 3.4 is of the same order as similar results in the literature, e.g. Suggala et al (2020) Corollary 4.3. We will make this clearer in the revised manuscript.

Minor comments:

Notational confusion: $\nabla_i, \Delta_i, \partial_i$: Our original motivation for using $\Delta_t$ stemmed from the connection to the discrete time step in the Euler discretization view of ResNets (Eq. 1). Recognizing the ambiguity, we have replaced $\Delta_t$ (the linear map within the residual block) with $A_t$ throughout the paper. We have also ensured that the text in Section 2 and Appendix B/C clearly specifies when partial derivatives or gradients are taken with respect to specific function arguments.

Line 111 $W_t$ shape: Thanks for spotting this. We will correct the dimensions in the revised manuscript.

Lines 204-207 'diagonal $\Delta_t$': We have improved the text to make it clearer that we only consider $\Delta_{\text{diag}}$ and $\Delta_{\text{scalar}}$ in the case when hidden size is equal to feature dimension size $D=p$. The dense case is free of this restriction, allowing for rectangular matrices.

Sandwiched Least Squares: We have modified the paper to properly address the 'sandwiched least squares' problems as Generalized Sylvester Equations.

Usage of 'amount of say': While "amount of say" is standard terminology used in the classical boosting literature, notably in AdaBoost [Freund & Schapire, 1997], we understand it may be unfamiliar to a broader audience. We have adopted your suggestion and replaced it with "line search step size" or similar phrasing for improved clarity.

---

We hope these responses and revisions address your concerns adequately.