Papers

Accurate energy demand forecasting is essential for the reliable operation and planning of modern sustainable energy systems. Spatial-temporal graph neural networks (STGNNs) have recently achieved strong performance in point forecasting by jointly modeling temporal dynamics and relational dependencies across interconnected energy nodes. However, in real-world energy systems, accurate point forecasts alone are insufficient, as operators also require reliable uncertainty estimates to support risk-aware decision-making, grid stability, and operational planning under uncertainty. Conformal prediction provides a principled and model-agnostic framework for uncertainty quantification under exchangeability assumptions, making it particularly attractive for safety-critical energy applications. However, existing conformal prediction approaches often fail to fully capture the complex spatial-temporal structure of energy systems. To address these limitations, we propose STOIC (Spatial-Temporal uncertainty quantificatiOn via In-Context learning), a conformal-inspired post-hoc residual calibration framework that integrates graph-based forecasting with the in-context learning capabilities of tabular foundation models. STOIC first generates point forecasts using an STGNN and subsequently reformulates spatial-temporal residuals into a tabular representation suitable for in-context learning. Leveraging a tabular foundation model, STOIC estimates feature-conditional residual quantiles without task-specific retraining of the calibration model, effectively capturing both sequential and relational dependencies. STOIC should be viewed as a conformal-inspired residual calibration method and does not by itself provide a finite-sample distribution-free coverage guarantee under arbitrary temporal dependence...

We study policy optimization for online episodic tabular Markov decision processes with unknown transition kernels, aiming for best-of-both-worlds guarantees together with data-dependent regret bounds. Recent work (Dann et al., 2023; Li et al., 2026) has shown that policy optimization can adapt to both adversarial and stochastic losses with first-order, second-order, and path-length bounds, but only under known transitions, leaving open whether such data-dependent guarantees are achievable by policy optimization when the transition kernel is unknown. We resolve this by developing a new algorithm based on optimistic follow-the-regularized-leader that attains these guarantees under unknown transitions. The key ingredient is a new design of optimistic $Q$-function estimators together with a data-dependent transition bonus that controls estimator bias through the loss-prediction error. Our analysis further identifies an unavoidable transition-dependent complexity term that captures the intrinsic cost of estimating the transition kernel. As a result, we obtain first-order, second-order, and path-length bounds with the transition-dependent complexity term while simultaneously achieving gap-dependent $\mathrm{polylog}(T)$ regret in the stochastic regime.

Nitrogen fertilizer management plays a central role in balancing agricultural productivity and environmental sustainability, yet identifying optimal application strategies remains difficult because treatment responses vary substantially across locations and many fertilizer choices are statistically indistinguishable near the optimum. This paper develops a hierarchical refinement procedure, built on sequential screening, for fertilizer recommendation in multi-site experiments that explicitly accounts for spatial heterogeneity while prioritizing parsimonious, decision-oriented selection. Rather than targeting a single estimated best treatment, the proposed method first conducts sequential screening at a higher aggregation level to eliminate clearly inferior fertilizer choices and then refines recommendations locally among the surviving candidates. We study the asymptotic properties of the proposed estimators and show that it provides screening-safety guaranteed recommendations. The efficacy of the new approach is investigated through a multi-state, multi-year corn nitrogen trial. The results show that no single fertilizer regime is uniformly optimal within a state; instead, each state is associated with multiple recommended choices, and the most common recommendation typically covers only about one-third to one-half of decision units, underscoring substantial within-state heterogeneity. Representative site-level comparisons further demonstrate that the proposed method often yields lower total nitrogen recommendations than state-level or hindsight benchmarks while maintaining competitive agronomic performance.

The Gaussian Kernel Robust Regression method (GKRReg) is a robust regression estimator that iteratively re-weights observations via a Gaussian kernel so that outliers and leverage points receive near-zero weight, with convergence of the estimation algorithm theoretically guaranteed. Despite a thorough study of estimation, the original work leaves open the problem of statistical inference for the regression coefficients. We fill this gap with three contributions. First, we formally establish that GKRReg belongs to the family of redescending M-estimators, providing the theoretical foundation for the inferential procedures that follow. Second, we derive a closed-form analytic sandwich variance estimator based on the theory of generalised M-estimators, corresponding to the HC0 class of heteroskedasticity-robust covariance matrices; we show that a finite-sample correction analogous to HC3 requires the weighted hat matrix of the converged IRWLS step, and identify this as a direction for future work. Third, we propose a pairs bootstrap that re-estimates the kernel width hyper-parameter gamma^2 on every replicate, capturing variability that the sandwich ignores. All procedures are implemented in the R package gkrreg, which also provides four estimators for gamma^2 and an automatic data-driven selection procedure, comprehensive diagnostic plots, and six real datasets from the robust regression literature. Applications to real data sets and comparison with traditional robust regression models highlight the potential of the GKRReg and the usability of the R package.

We consider the calibration of probability forecasts. Several notions of calibration exist when the forecaster issues a single forecast for each of the observations that is to be predicted. We extend one of these notions, auto-calibration, to the common situation in which the forecaster issues a sequence of forecasts for each observation, repeatedly updating their prediction as they receive additional information. For observations that sit in any Borel space, we show that auto-calibration is equivalent to a certain sequence of random probability measures satisfying the martingale property, and we propose a simple, statistical approach to testing this property. This provides, for the first time, a way of testing the calibration of such sequences of probability forecasts.

Conformal prediction and its variants, including the split conformal prediction, provide a distribution-free framework for uncertainty quantification by constructing prediction intervals or sets with finite-sample coverage guarantees. The statistical efficiency of these intervals depends critically on how the data are split into training and calibration samples. Despite its practical importance, a principled characterization of the training-calibration split that minimizes prediction interval length while maintaining coverage has remained largely unresolved. In this paper, we develop a theoretical framework for optimal data splitting in split conformal prediction. We first analyze the problem in a general setting and derive analytical characterizations of the length-optimal split ratio under both symmetric and asymmetric regimes. We then show how the general results specialize to several commonly used regression settings, including linear regression, nonparametric regression, and neural networks, thereby demonstrating the scope of the framework. We also describe a data-based method for selecting the optimal proportion. Our analysis clarifies how model-related features govern the optimal allocation of samples between training and calibration and provides principled guidance for constructing shorter prediction intervals. Experiments on both synthetic and real-world datasets demonstrate the applicability of the proposed methodology across a variety of practical scenarios.

Classical minimax lower bounds for testing are typically derived for fixed error probabilities, while high-confidence results often impose a common failure probability. We study prescribed-error testing, in which the level and the target type II error may be small and of different orders. Standard prior-based reductions generally aggregate the two errors into a single quantity and therefore do not capture their distinct roles. We develop a general lower-bound technique based on a binary reduction that preserves the separate roles of the two error targets. The reduction yields two directed Kullback-Leibler information requirements, corresponding respectively to the level and the target type II error. When both directed mixture divergences can be controlled, they combine into a binary Jeffreys divergence, leading to the logarithmic dependence on the level and the target type II error. Applying the framework to Gaussian sequence testing, multinomial uniformity testing, and continuous uniformity testing over Hölder balls, we obtain lower bounds that match corresponding high-confidence upper bounds and hence establish prescribed-error minimax rates sharp up to constant factors.

The Self-Improving Alignment (SAIL) algorithm addresses distribution shift by reducing a bilevel formulation of the problem to an efficient, single-level method. Empirically, SAIL has demonstrated strong performance on this task. However, a formal analysis of its convergence properties has been lacking. We identify a key theoretical challenge: the standard SAIL objective function is not guaranteed to be strongly concave due to unfavorable properties of its Hessian. To address this limitation, we propose a regularized objective, SAIL-RevKL, which incorporates a reverse Kullback-Leibler (KL) divergence penalty to improve the optimization landscape. Our central theoretical contribution is to prove that this regularized objective satisfies the Polyak-Lojasiewicz (PL) condition within a bounded parameter space. We establish global convergence guarantees, achieving a near-linear sample complexity. We further validate the effectiveness and stability of SAIL-RevKL through empirical evaluations, demonstrating that it outperforms the vanilla SAIL on both MuJoCo benchmarks and LLM alignment tasks.

Functional principal component analysis has traditionally focused on processes indexed by Euclidean domains, although functional observations increasingly arise on curved domains. We develop FPCA for real-valued processes observed with noise at discrete locations on a general Riemannian manifold. We estimate the mean and covariance functions using intrinsic kernel smoothing based on geodesic distance and a normal-coordinate volume correction, with a general weighting scheme accommodating unequal numbers of observations across subjects. We establish $L^2$ and uniform convergence rates for these estimators and derive corresponding error bounds for the estimated eigenvalues and eigenfunctions. These results quantify how manifold dimension, sampling frequency and weighting affect estimation accuracy and characterize the phase transition between sparse and dense sampling regimes. Simulation studies on the circle and sphere, together with real-data applications, illustrate the finite-sample performance of the proposed methodology.

We investigate the contextual slate bandit problem with generalized linear rewards under limited adaptivity. At each round, the learner is presented with $N$ sets of items, where each item is represented by a $d$-dimensional feature vector. The learner then constructs a slate by selecting one item per set; the resulting slate yields a scalar reward sampled from a Generalized Linear Model (GLM). We propose algorithms under two limited-adaptivity settings: (a) Batched and (b) Rarely-Switching. For the batched setting, we introduce B-SlateGLinCB, which partitions the time horizon into $\mathcal{O}(\log\log T)$ batches such that each batch's policy relies only on data from previous batches. For the rarely-switching setting, we propose RS-SlateGLinCB, which adaptively performs only $\mathcal{O}(Nd\log T)$ parameter updates. Under a diversity assumption on the item sequences, we prove that B-SlateGLinCB and RS-SlateGLinCB achieve regret bounds of $\mathcal{O}(Nd^{3/2}\sqrt{T})$ and $\mathcal{O}(Nd\sqrt{T})$, respectively. Notably, both bounds are independent of the non-linearity parameter $\kappa$ that is typically found to scale the regret of GLM bandit algorithms. Our algorithms are computationally efficient, requiring only $\text{poly}(N)$ time per round despite $2^{\Omega(N)}$ possible slates. Simulations show our algorithms outperform existing baselines with limited adaptivity and remain competitive with Slate-GLM-OFU, a fully adaptive state-of-the-art algorithm. Notably, a slightly modified B-SlateGLinCB empirically matches this baseline. Finally, we demonstrate strong performance in a practical in-context example selection task for language models.

We introduce a geometric formulation of statistical feature learning for supervised regression. Feature learning is defined through a base–fiber decomposition: the base is the feature-side geometry produced by training, and the fiber is the learned feature space where estimation is performed. We prove this property for spherical mean-field Langevin dynamics, viewed as the Wasserstein gradient flow of a negative entropy-regularized empirical risk. In Gaussian multi-index models, the low-temperature stationary distribution concentrates near the hidden indices, forms a multi-spike structure, and yields parameter recovery with high probability, even though negative entropy regularization penalizes concentration. This concentration has a sharp transition at temperature $\lambda\asymp 1$. In Gaussian single-index models, the stationary measure satisfies a concentration property, with parity determining whether it lives on $S_2^{d-1}$ or $\mathbb{RP}^{d-1}$. The induced learned feature space aligns the regression signal and yields rates $d/N$ and $Md/N$, up to logarithmic factors.

Demographic corrections are routinely performed in many disciplines, including psychology. Yet, there are ongoing debates about whether these corrections are appropriate and improve classification accuracy. Here, we focus on cognitive screening tests, and show that common demographic corrections, like the z-score correction, can be detrimental for classification in some settings. Formally, we present sufficient conditions ensuring that raw scores outperform the demographically-corrected ones, and give a substantive interpretation of this result. We also investigate the claim that using demographically-corrected scores results in more fair decisions compared to using raw scores. We apply our results to the Mini-Mental State Examination in the OASIS-3 dataset.

The increasing availability of diverse data sources has motivated great interest in data integration for improving regression efficiency. Existing data integration methods primarily focus on integrating nonprobability samples and typically assume that the integrated data sources represent the same target population. While this assumption is often difficult to justify for nonprobability samples, it is naturally satisfied when integrating probability-based surveys designed to represent a common target population. Such surveys are important research data sources because they provide representative samples and collect rich information on diverse variables, making them well suited to data integration. However, existing integration methods do not accommodate complex sampling designs. We propose model-assisted calibration methods to improve regression efficiency by integrating multiple probability-based survey samples. The proposed framework accommodates settings in which either individual-level data or only summary statistics are available from external surveys while preserving valid finite-population inference without requiring correct specification of the outcome model. We establish the design consistency of the proposed estimators and develop Taylor linearization variance estimators accounting for the complex sampling designs of both surveys. Simulation studies and an application integrating National Health and Nutrition Examination Survey and National Health Interview Survey demonstrate substantial efficiency gains while maintaining valid finite-population inference.

We derive sufficient and almost necessary conditions for large sample and catalytic majorization between finite statistical experiments over standard Borel sample spaces. This work generalizes previous results, on one hand, in the bivariate case and, on the other hand, in the multivariate discrete (or, rather, finite) case, i.e., matrix majorization. We derive multivariate generalizations of the bivariate Renyi relative entropies and show that inequalities involving these multivariate Renyi divergences characterize large-sample and catalytic majorization of finite statistical experiments. As our methods are real-algebraic in nature, this work demonstrates that large deviation techniques are not the only option available to derive conditions for large sample majorization even in the case of more general sample spaces of the experiments. We also show that all general multivariate divergences, i.e., multivariate extensive and monotone maps of finite statistical experiments, can be expressed through barycentres over the set of multivariate Renyi divergences. We also show that we may characterize the optimal conversion rate of a statistical experiment into another using the multivariate Renyi divergences.

Quantile regression aims to estimate the conditional quantiles of a response variable from observed data. In a Bayesian setting, Gaussian process quantile regression provides uncertainty quantification but faces significant computational challenges due to the nonconjugacy of the asymmetric Laplace likelihood and the cost of posterior inference. We develop a sparse Gaussian process framework in which the quantile function is represented through a reduced set of inducing variables and posterior inference is performed using a Laplace approximation. A decomposition of the predictive uncertainty into conditional-prior and posterior-induced variance components is then exploited to drive two complementary adaptive mechanisms: inducing-input infilling and data acquisition. These mechanisms are combined within a sequential algorithm that allocates computational effort toward the dominant source of predictive uncertainty and adaptively controls model complexity. Numerical experiments on benchmark problems demonstrate the accuracy of the Laplace approximation, the benefits of variance-based inducing-input placement, and the effectiveness of the proposed sequential enrichment strategy compared with predefined data-acquisition strategies. The method is finally applied to ship seakeeping under forecast uncertainty, where it outperforms a replication-based counterpart from the literature under limited simulation budgets.

Functional data analysis (FDA) provides statistical methods for analyzing samples of time-continuous stochastic processes. Measurements often arise in the form of sensor data for a key scientific variable. The practical problem of irregular sensor disruptions has fostered interest in analyzing partially observed random functions. Specifically, this paper is motivated by a time series of intermittently missing pollution data with dependence along pollution paths and missingness patterns. To allow statistical analysis, we develop the first inference methods for dependent, partially observed functional time series. Existing methods were not appropriate for this task, because they heavily rely on the independence of the data functions. Mathematically, we model data on the space of bounded functions equipped with the supremum norm. This allows simultaneous inference across the entire functional domain, including simultaneous confidence bands – something existing Hilbert-space-based methods cannot provide. To study non-stationary trends along the time series, we extend state-of-the-art multiscale inference methods (originally developed for scalar data) to partially observed functions. The key application of the latter methods is testing for excessive pollution levels in inner cities. Our approach combines state-of-the-art Gaussian approximations with stochastic process theory. Interestingly, it also improves existing results for fully observed functional time series by avoiding a functional CLT.

The classical $k$-means clustering, based on distances computed from all data features, cannot be directly applied to incomplete data with missing values. A natural extension of $k$-means to missing data is to involve only the observed positions in clustering, which is equivalent to imputing missing values by corresponding cluster means. However, for data missing not at random (MNAR), since missingness is related to data values, such a mean-imputation-based method may lead to the distortion of estimated cluster centers, resulting in a poor clustering result. Since MNAR mechanisms are very common in reality, it is necessary to improve the performance of $k$-means-based clustering methods for such data. In this paper, we focus on a magnitude-decaying MNAR scenario where data is more likely to be missing at positions with smaller absolute values, and we propose a novel $k$-means clustering method based on the constraint of the size of imputation values, which enjoys a good mathematical interpretation. Moreover, we establish the statistical consistency of the estimated cluster centers of the proposed method to the true cluster centers of fully observed data, and solve the optimization of the proposed loss function via an alternative minimization algorithm. Simulation experiments verify the effect of the proposed method in improving clustering results and reducing the bias of estimated cluster centers. Applications to real-world missing data further show the utility of the proposed method.

We study the task of learning the structure of a $d$-sparse Gaussian graphical model on $n$ variables from a single trajectory of Glauber dynamics. Beyond algorithmic considerations, many applications present temporally correlated observations rather than i.i.d.\ samples. In the classical i.i.d.\ setting, under comparably general sparsity and minimum edge-strength assumptions, sublinear-in-$n$ sample guarantees are known, but achieving them in polynomial-time remains open. Motivated in part by this gap, we give a polynomial-time algorithm that recovers the conditional-independence graph from a single Glauber trajectory, with a trajectory-length guarantee that does not depend on the mixing time. Technically, our algorithm has three components. First, we estimate the conditional variances and rescale the trajectory to reduce to the unit-diagonal case, without changing the underlying graph. Second, we design a local edge test that extracts adjacency information from short update windows by isolating pairwise influence. Third, we aggregate these local statistics using a robust median-based estimator, and prove accuracy despite temporal dependence arising from a single trajectory.

Modern experiments, including evaluations of AI-enabled services and platform interventions, often depart from independent and identically distributed (i.i.d.) sampling because assignments may be adaptive, balanced across covariates, or subject to rollout constraints such as exposure, fairness, and budget limits. This paper studies the efficiency benchmark for estimating causal targets in such sequential experiments. We show that every non-anticipating design induces an average propensity score, and we establish a semiparametric lower bound: for regular locally unbiased estimators, attainable precision is bounded by the i.i.d. efficiency benchmark evaluated at this induced score. The average propensity score thereby serves as a common benchmark and design target, allowing sequential experimental design to be viewed as choosing or learning an efficient allocation rule, with operational constraints entering through the admissible set when present. We then develop implementable batched adaptive designs that approach this benchmark through two complementary mechanisms. The first uses regression adjustment based on efficient influence functions; for general smooth estimands it attains the benchmark under standard nuisance-rate conditions, while for linear functionals of outcome means it achieves a sharp second-order rate. The second uses adaptive covariate balancing to attain the same benchmark through the assignment mechanism, enabling simple moment-based estimation. Both routes require only a small number of policy updates, making them compatible with delayed feedback and easier to monitor in operational deployments. Numerical experiments and an empirical study of AI medical-assistant evaluation demonstrate the practical efficiency gains, including in multi-treatment settings. Overall, the paper provides a unified framework for characterizing and designing efficient sequential experiments.

Predicting biomolecular properties from limited labeled data is a central bottleneck in protein engineering and small-molecule design. As strong pretrained encoders now supply rich fixed-length representations, the difficulty has shifted from representation learning to building a data-efficient predictor for the few-shot regime. Tabular foundation models such as TabPFN and TabICL are unlikely candidates for this role: they are in-context learners pretrained on synthetic tables drawn from random causal graphs, a generative prior with no obvious correspondence to the processes that produce protein sequences or molecular graphs. That this tabular, causal inductive bias should transfer to biomolecular data at all is counter-intuitive, yet we find it does. Treating each method as a predictor-representation pair, we evaluate across two domains. We find that on protein fitness regression tasks these in-context learning models coupled with ESM Cambrian representations achieve or exceed state-of-the-art results on ProteinGym, and outperform task-specific supervised regressors on a diverse esterase catalytic activity dataset. For small-molecule classification with ECFP/RDKit descriptors, no single predictor-representation pairing dominates across TDC ADMET, MoleculeNet, FS-Mol, and DrugOOD, but they are competitive with the existing task-specific state-of-the-art. Crucially, on both protein and small-molecule few-shot tasks, these predictor-representation pairs offer strong performance. We conclude that tabular foundation models can be strong biomolecular predictors, but only when coupled with expressive representations.

Broken adaptive ridge (BAR) penalty approximates $L_0$-regularization through iterative reweighting of L2 penalties. This penalty enjoys both the oracle property and the grouping effect for highly correlated covariates, making it particularly attractive for penalized regression with complex dependence among predictors. In this paper, we develop a BAR-penalized linear rank regression method for the semiparametric accelerated failure time model with right-censored data. Computational tractability is achieved by applying induced smoothing to the nonsmooth Gehan-type rank estimating function, yielding a more stable framework for estimation and inference. For scalable penalization, we develop a cyclic coordinate descent algorithm that minimizes the penalized objective function, and estimates the regression coefficients in a coordinate-wise manner. We further extend the proposed method to more complex survival endpoints, such as multivariate partly interval-censored (PIC) data. Under mild conditions, the proposed estimator satisfies both the oracle property and the grouping effect, and the variance estimator of the informative coefficients can be derived in analytic form. Numerical studies using synthetic data compare our approach to several well-known penalties, and demonstrate its superior selection accuracy and estimation efficiency across various scenarios. Furthermore, applications to right-censored outcomes from primary biliary cirrhosis, and correlated PIC outcomes from colorectal cancer further illustrate the practical utility of the proposed method. The R package aftPenCDA for implementing the method is available on R CRAN.

Gaussian process inference is often limited by cubic computational costs, a challenge that becomes more pronounced in spatio-temporal settings where posterior inference is required over dense grids. While state-space SPDE formulations enable linear complexity in time, exact inference remains cubic in space and deteriorates further when observation locations are disjoint from the prediction locations, which inflates the number of considered spatial points. To address this, we propose the Vanilla-SPDE Exchange, which exploits an equivalence between the standard and SPDE formulations of GP inference to construct a hybrid scheme with improved computational cost. We demonstrate these gains through complexity analysis and numerical experiments.

We study high-frequency inference for an ergodic regime-switching jump diffusion with locally stable infinite-variation jumps. The state and the finite-state regime are observed on the same grid. The mesh satisfies \(h_n\to0\), while the horizon satisfies \(T_n=nh_n\to\infty\). We estimate structural drift and diffusion parameters in the presence of unknown regime-specific Lévy measures. For a common activity index \(1<\beta<8/5\), we remove two leading jump biases from each regime's truncated variation. We then combine the diffusion estimator with estimated signed tail means in an adaptive drift contrast. We prove joint asymptotic normality at rates \(\sqrt{T_n}\) and \(\sqrt n\), together with feasible studentization. The growing horizon requires uniform control of the debiasing remainders and random denominators over increasing regime exposures. It also produces an explicit Lévy-measure term in the drift covariance. A fixed-grid median stabilizes the two-step correction without changing its first-order limit. We study finite-sample behavior by simulation and apply the method to 127.5 hours of cattle sensor data. In the application, jaw-activity labels are recorded by a device separate from the sensor that measures the state process.

Ranked choice voting (RCV) is a popular alternative voting method in which voters are asked to list their favored candidates in preference order, rather than vote for a single candidate. When these ballots are tabulated, candidates are successively eliminated, and their votes are reallocated to each voter's next-preferred choice. The process continues until a candidate commands a majority of the active ballots and is declared the winner. As RCV gains wider adoption, the method poses novel challenges for pollsters. Unlike plurality elections, the event that a candidate wins cannot be expressed in terms of a single population parameter. Hence, the basic concept of a margin-of-error is not straightforward to define. Moreover, a candidate's ability to win may depend on both their support across the ballot and the order in which other candidates are eliminated. Existing measures of sampling uncertainty for polls of RCV elections do not clearly quantify these path-dependent outcomes. Here, we propose a simple, Bayesian framework to quantify uncertainty in polls of RCV elections. We cast the problem as one of estimating win probabilities for each leading candidate, and leverage a simple conjugacy relationship to estimate these probabilities conditional on the poll results. We include applied analyses involving two prominent ranked choice voting elections: the 2021 New York City Democratic mayoral primary, in which Eric Adams narrowly defeated Kathryn Garcia in the final round; and the 2022 special election to Alaska's U.S. House seat, in which Mary Peltola was elected despite not being a Condorcet winner. Using the cast vote records from both elections, we demonstrate some challenges of traditional frequentist uncertainty quantification in RCV polls. We also demonstrate the utility of our approach using a poll of the NYC primary obtained from the polling firm Data for Progress.

In two-sided marketplaces with heterogeneous products, it is important to understand the causal relationship between additional supply and marketplace outcomes, such as the total quantity transacted or transaction value in the marketplace. This paper studies a causal machine learning approach to estimating this relationship across product segments. We use the Airbnb marketplace as an example, focusing on the impact of additional listing supply on total bookings, but the methodology applies to other two-sided marketplaces. Our approach combines double/debiased machine learning with a hierarchical Bayesian framework that leverages pre-existing knowledge as priors. We construct tractable and informative features for the model by leveraging measures of product segment similarity from the geospatial literature. We find that such a model provides plausible estimates of the marketplace returns to additional supply and strong out of sample performance.

Recent work has shown that well-optimized individual decision trees can match complex black box models in some settings, primarily in noisy domains. For the remaining settings, however, complex ensembled compositions of trees often achieve higher accuracy at the cost of interpretability, leaving practitioners with difficult modeling decisions along an accuracy-interpretability tradeoff. Ideally, we would like to classify as much of the data as possible with one or a small number of trees, achieving interpretability for most samples while maintaining state-of-the-art accuracy. We introduce Multistage Defer Trees: a sequence of sparse decision trees that each make predictions for most samples, while deferring a small proportion to the next tree in the sequence or, ultimately, to a black box. We demonstrate that we can train this model class to match the performance of complex tree-based ensembles while routing most samples through only one or a small number of sparse decision trees. We discuss a range of techniques for training these models while maintaining simplicity. Our method expands the accuracy–interpretability frontier in settings where single-tree methods remain insufficient, demonstrating that even when complex models are necessary, they need not be fully opaque.

In biological data from allometry studies, the largest eigenvalue is typically dominant, and the gaps between minor eigenvalues are often narrow. Such proximity among small minor eigenvalues can lead to instability in statistics based on their corresponding eigenvectors. This study derives the asymptotic normality of the hybrid principal component analysis estimator of the leading principal eigenvector in the multivariate allometric regression model and proposes a test based on a geometric statistic for the parallelism between the regression direction and the principal component direction that avoids this instability. Using the hybrid principal component analysis framework, we analyze the well-known painted turtle carapace data and confirm previously reported results on the allometric extension relationship between female and male turtles.

Multicollinearity is a long lasting challenge in observational causal inference, especially in regressions – highly correlated independent variables make it hard to isolate their individual impacts on outcomes of interest. While common solutions such as shrinkage estimators and principal component regressions are helpful in prediction problems, a crucial limitation hinders their applicability to causal inference problems – they cannot provide the original causal relationships. To fill the gap, we present an innovative and intuitive solution, by employing hierarchical clustering to aggregate data in a way that effectively alleviates collinearity. This method is generally applicable to causal problems featuring multicollinearity. We use a marketing application to demonstrate how and why it works. Expenditures on different advertising channels often exhibit correlations, making it exceedingly difficult to separately measure their impact. Many previous studies proposed to leverage granular cross-sectional data for better identification but, to our knowledge, none explicitly addressed multicollinearity, which undermines causal identification even with granular data. We propose to hierarchically cluster geographic units based on marketing spend correlation to reduce collinearity, and to implement a Bayesian Marketing Mix Model with cluster-level data. Such clustering happens in two steps – we first normalize and demean geo-level data to establish a common scale and to eliminate the common trends; we then calculate pairwise distance to summarize marketing spend correlation between geos and cluster the ones with moderate to strong correlation. Both descriptive evidence and regression analysis affirm that such hierarchical clustering effectively mitigates collinearity and facilitates the separate identification of the impact of different marketing channels.

The study of the Earth's behavior has greatly benefited from the widespread deployment of Global Navigation Satellite Systems (GNSS), enabling large-scale monitoring of crustal deformation and long-term geophysical trends. In this work, we focus on the North American region, where complex tectonic activity, particularly along the western margin, requires methods capable of processing and analyzing large collections of GNSS time series distributed across extensive spatial domains. Analyzing the full GNSS network provides a coherent view of deformation across multiple scales, allowing detection of long-wavelength signals, subtle intraplate strain and regionally consistent velocity fields that are difficult to capture through local or subsampled analyses. Despite the availability of such data, existing methodologies remain computationally prohibitive for large-scale analyses across this region (or others). To address this limitation, we introduce a highly scalable framework (implemented in open-source software) that enables inference on crustal deformation across the full North American GNSS network using standard computational resources. The proposed method achieves substantial computational gains while maintaining inferential performance comparable to existing approaches and confirms existing tectonic trends, thereby supporting efficient large-scale monitoring and contributing to ongoing efforts toward "Open Science".

Resampling-based simultaneous confidence bands for cumulative hazard functions often undercover in finite samples with right censoring. We study two aspects of the construction that can contribute to this gap, the resampling scheme and the calibration statistic, and propose a procedure that intervenes on both. The exchangeable bootstrap reweights the numerator and the denominator of the Nelson-Aalen ratio, preserving its ratio structure. The box-calibrated discrepancy constructs lower and upper step envelopes from adjacent values of the original and resampled Nelson-Aalen estimators and measures the resulting vertical discrepancy. We establish conditional weak convergence of the exchangeable bootstrap, prove that box calibration is first-order asymptotically equivalent to grid calibration, and show that the resulting band attains nominal coverage asymptotically. The box correction uses the same bootstrap paths and event-time grid as grid calibration; after each bootstrap path is formed, it requires only an additional linear pass over the event-time grid and therefore has negligible computational overhead. In simulations across a range of hazard shapes and censoring levels, the exchangeable bootstrap with box calibration is, in most configurations, closest to nominal coverage among the methods considered. A notable consequence is a ranking reversal: the ratio-preserving exchangeable bootstrap has the lowest coverage under grid calibration, yet is usually closest to the nominal level after box calibration. A melanoma data example illustrates the practical effect on the cumulative hazard bands. The proposed procedure operates on the original cumulative-hazard scale, requires no variance-stabilizing transformation, and permits inference from time zero.