PRTree: Smooth Probabilistic Regression Trees

• PRTree Package is an R tool for fitting smooth regression trees that replace hard splits with probabilistic region memberships.
• It efficiently handles incomplete predictor data by employing three strategies—uniform probability, partial observation, and dimension-reduced smoothing—without external imputation.
• The package integrates R, FORTRAN, and C to deliver rapid model fitting with SVD-based least squares and supports multiple kernels like Gaussian and Student’s t.

Searching arXiv for the cited PRTree-related papers to ground the article in the latest relevant literature. arXiv search: (Prass et al., 4 Oct 2025) PRTree is an R package for fitting Probabilistic Regression Trees (PRTrees), a smooth generalization of classical regression trees in which hard, axis-aligned region membership is replaced by probabilistic region memberships. In the package formulation, PRTrees are extended to incomplete covariate data through three implemented strategies—uniform probability, partial observation, and dimension-reduced smoothing—so that missing predictor values can be handled internally rather than through external imputation. The package is distributed on CRAN and is implemented with R + FORTRAN + C (Prass et al., 4 Oct 2025).

1. Mathematical formulation

The package is built around the regression setting

$$Y = f(X) + \varepsilon,\qquad E(\varepsilon)=0,\quad E(\varepsilon^2)<\infty.$$

Its starting point is the contrast between CART and PRTrees. CART approximates the regression function by partitioning the feature space into disjoint regions $R_1,\dots,R_M$ and predicting

$$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$

PRTrees replace the indicator $I(X\in R_m)$ by a smooth probabilistic association function $\Psi$, yielding

$$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$

where $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$, $R_m$ are tree regions, $\gamma_m$ are leaf or region coefficients, and $\sigma \in \mathbb{R}^p_+$ controls smoothness. The association function is defined by convolution over each region: $R_1,\dots,R_M$0

Two structural properties are emphasized. First, $R_1,\dots,R_M$1 acts as the probability-like weight that observation $R_1,\dots,R_M$2 belongs to region $R_1,\dots,R_M$3. Second, the leaf weights satisfy

$R_1,\dots,R_M$4

This preserves the tree interpretation while replacing piecewise-constant prediction by continuous prediction. The paper states that CART is recovered as a limiting or hard-split special case, so the package is best understood as a smooth tree model rather than a departure from tree structure itself (Prass et al., 4 Oct 2025).

2. Missing-data recursion and the three implemented methods

The extension to incomplete data is formulated under MCAR. For a subset $R_1,\dots,R_M$5, $R_1,\dots,R_M$6 denotes the subvector of observed coordinates in $R_1,\dots,R_M$7. For a rectangular region $R_1,\dots,R_M$8,

$R_1,\dots,R_M$9

The original smoothing function induces a probability measure

$$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$0

The paper highlights two useful properties: marginal compatibility for unbounded coordinates, and probability conservation under splitting,

$$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$1

The package defines a proxy function $$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$2 and recursively constructs a missing-data-aware probability: $$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$3 At the root,

$$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$4

A central rule is that if the splitting variable itself is missing, then

$$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$5

so the observation is split equally between the two branches.

The three implemented strategies are controlled by fill_type:

fill_type	Definition of $$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$6	Stated effect
0	$$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$7 if $$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$8 is fully observed; otherwise $$f_{\mathrm{CART}}(X)=\sum_{m=1}^M c_m\, I(X\in R_m).$$9	Uniform weight under incompleteness
1	$I(X\in R_m)$0 for $I(X\in R_m)$1	Hard compatibility on observed coordinates
2	$I(X\in R_m)$2 for $I(X\in R_m)$3	Smooth probabilistic routing on observed coordinates

The uniform probability method (fill_type = 0) is described as the most conservative and least informative method. The partial observation approach (fill_type = 1) uses observed values to eliminate incompatible regions, but applies smoothing only when the observation is fully complete. The dimension-reduced smoothing technique (fill_type = 2) projects both the observation and the region onto the observed coordinates and applies the original PRTree smoothing there. The paper characterizes this third strategy as the most faithful extension of the original PRTrees to incomplete data, and the simulations identify it as the best-performing missing-data strategy among the three (Prass et al., 4 Oct 2025).

A computationally important product representation is also given. If $I(X\in R_m)$4 is the set of internal nodes along the path from the root to $I(X\in R_m)$5, then

$I(X\in R_m)$6

with

$I(X\in R_m)$7

The paper notes the special case that if $I(X\in R_m)$8 and $I(X\in R_m)$9 is missing, then

$\Psi$0

where $\Psi$1 is the depth of the node.

3. Estimation, split search, and stopping rules

Given training data $\Psi$2, parameter estimation is posed as empirical risk minimization with squared loss: $\Psi$3 For fixed regions and $\Psi$4, the coefficient vector is estimated by ordinary least squares,

$\Psi$5

when $\Psi$6 is nonsingular. The implementation does not explicitly form the inverse; instead it solves the least-squares problem using LAPACK DGELSD via an SVD-based divide-and-conquer routine, which the paper presents as numerically stable even when $\Psi$7 is singular or nearly singular.

Tree growth is greedy. For each current leaf or node, the algorithm searches over candidate split variables $\Psi$8 and thresholds $\Psi$9, and the selected split minimizes the global MSE after the tree is updated: $$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$0 The paper emphasizes a two-stage splitting algorithm. In the first stage, candidate splits are scored using the proxy criterion

$$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$1

Candidates may be selected by node or globally across nodes. During this stage, missing values are temporarily assigned using proxy_crit, with options "mean", "var", and "both".

In the second stage, each retained candidate split is fully evaluated by updating the tree, recomputing $$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$2 and $$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$3, computing global MSE, and selecting the split with the largest MSE reduction. A split is accepted only if the improvement exceeds the complexity threshold cp. The package also uses max_depth, max_terminal_nodes, n_min, and the pair perc_x and p_min as stopping rules. The criterion associated with cp is stated as

$$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$4

The smoothness vector $$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$5 is selected by grid search: the package computes sample standard deviations $$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$6, builds a default grid of candidate vectors by scaling $$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$7, evaluates each candidate on a validation set, and selects the $$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$8 giving the lowest MSE. By default, grid_size = 8 and the multipliers range from 0.25 to 2.00 (Prass et al., 4 Oct 2025).

4. Software architecture and package interface

The implementation is divided across three layers: R for user interface and orchestration, FORTRAN for tree construction and iterative search, and C for optimized probability computations. The package supports several kernels through the probability measure $$f_{\mathrm{PR}}(X;\Theta)=\sum_{m=1}^{M}\gamma_m\,\Psi(X;R_m,\sigma),$$9: Gaussian (norm), Log-normal (lnorm), Student’s t (t), and Gamma (gamma). The paper associates these with parameters such as sdlog, df, and shape (Prass et al., 4 Oct 2025).

The principal control constructor is pr_tree_control, which creates and validates the control object and returns an object of class prtree.control. The paper lists the main configurable arguments as cp, max_depth, n_min, fill_type, proxy_crit, n_candidates, by_node, dist, dist_pars, max_terminal_nodes, and others. This control structure exposes the main modeling choices: regularization, depth and size constraints, missing-data routing, candidate generation, and kernel family.

The fitting function is pr_tree. Its inputs are y, a numeric response vector; X, a numeric matrix or data frame of predictors; control, usually from pr_tree_control; and optional direct control arguments. The described preprocessing pipeline verifies that y has no missing values, processes training and validation indices, can perform a stratified split if idx_train is not supplied, and generates sigma_grid automatically if none is given. The returned object of class prtree includes yhat, P, gamma, sigma, training and validation MSE, nodes_matrix_info, and regions.

Prediction is handled by predict.prtree. It computes the probability matrix $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$0 for new observations using the stored tree and applies the same distribution, $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$1, and fill_type used in fitting. With complete = FALSE, it returns predictions only; with complete = TRUE, it also returns the probability matrix. This makes the package suitable not only for point prediction but also for direct inspection of the soft region memberships that define the model’s interpretability.

5. Empirical behavior under MCAR

The simulation study in the paper uses the data-generating process

$\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$2

with $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$3, $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$4, and $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$5. Each replication contains $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$6 observations, with 1000 for training and 200 for testing. Missingness is MCAR at levels

$\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$7

For $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$8, $\Theta = (\{R_m\}_{m=1}^M,\gamma,\sigma)$9 training observations are selected, and each selected observation receives one of three patterns with equal probability: both $R_m$0 and $R_m$1 missing, only $R_m$2 missing, or only $R_m$3 missing. The same procedure is applied to the test data.

The comparison is against CART via rpart with default settings, and against PRTree with fill_type = 0,1,2. The PRTree configuration reported for simulation is grid_size = 8, max_terminal_nodes = 50, max_depth = 49, default cp = 0.01, n_min = 5, perc_x = 0.1, p_min = 0.05, proxy_crit = 3, by_node = FALSE, n_candidates = 3, and dist = "norm", with the first 800 training observations used to fit and the last 200 to select $R_m$4. Performance is measured using RMSE against both the noisy observed response $R_m$5 and the true regression function $R_m$6, on both training and test sets (Prass et al., 4 Oct 2025).

Several findings are stated explicitly. CART trees become slightly smaller as missingness increases. PRTree trees are generally smaller or comparable, and shrink much more under severe missingness. At $R_m$7, PRTree produces very parsimonious trees, with median terminal nodes around 4–6 depending on fill_type, versus about 7 for CART. Under complete data, PRTree has substantially lower true RMSE than CART: roughly 0.57–0.58 for PRTree versus 1.13–1.22 for CART. The paper also states that PRTree shows much better alignment between train and test error, indicating less overfitting. As missingness increases, all methods degrade, but CART deteriorates the most. Under high missingness ($R_m$8), PRTree remains clearly better, and fill_type = 2 is the best performer.

These results position the package as a smooth tree method whose empirical advantages are not restricted to incomplete-data settings. In the reported experiments, it is more accurate, more robust to missingness, more resistant to overfitting, and more parsimonious than CART (Prass et al., 4 Oct 2025).

6. Scope, limitations, and nomenclature

The intended use cases described for the package are regression problems with incomplete covariates, settings where interpretability matters, smooth nonlinear regression surfaces, applications where standard CART is too coarse, and workflows that benefit from an integrated R tool. The package’s practical advantages are stated as direct handling of missing predictors, preservation of the interpretability of trees, smooth and continuous predictions, support for multiple kernels and missing-data strategies, efficient implementation through FORTRAN, C, and LAPACK SVD routines, and return of rich tree structure and probability information for interpretation (Prass et al., 4 Oct 2025).

The limitations are equally explicit. The theory and simulations assume MCAR; MAR and MNAR are not handled theoretically. The tree remains a greedy model and is therefore not globally optimal. The choice of kernel and $R_m$9 may matter. The paper also notes that extending consistency theory to missing-data handling remains future work. The three missing-data strategies define a spectrum: fill_type = 0 may discard useful information, fill_type = 1 is hard and less smooth, and fill_type = 2 is the most effective method in the reported simulations.

A common source of confusion is that PRTree is not a unique acronym across the literature. In neural network interpretation, the name is used for the rank projection tree framework, which generates multiscale network interpretations and prioritized collections of subsets of inputs rather than only singleton rankings (Warrell et al., 2018). In Bayesian symbolic regression, PRTree or pRTree refers to a framework based on probabilistic regular tree priors, probabilistic Regular Tree Expressions, and probabilistic tree automata (Schneider et al., 2023). A plausible implication is that bibliographic searches for “PRTree” require explicit disambiguation between Probabilistic Regression Trees, rank projection trees, and probabilistic regular tree priors.

Within that broader naming landscape, the package discussed here is specifically the CRAN-distributed R implementation of Probabilistic Regression Trees with native missing-data handling. Its defining contribution is not a generic tree API, but a smooth tree estimator in which probabilistic region memberships, kernel-based routing, and missing-data recursion are integral parts of the model rather than post hoc modifications (Prass et al., 4 Oct 2025).