Papers

Quasi-Monte Carlo methods allow computer experiments to be run with far fewer simulations than crude Monte Carlo. A digital net in base two, such as Sobol', is however only balanced when the number of samples is a power of two, and Latin Hypercube Sampling (LHS) accepts any sample size but only controls the one-dimensional margins. This work proposes a space- filling design defined for any sample size, referred to as LAT for Latinized aperiodic tiling. The unit hypercube is cut recursively across its longest edge following the golden section into N cells of equal volume. One point is then placed in each cell, and the margins are made Latin while every point stays inside its own cell. The construction only uses integer splits and costs O(dNlog N). The recursion is shown to follow the Fibonacci word and its aperiodicity is analysed. LAT is assessed with four L2-discrepancies and with the integration error on analytical functions and engineering emulators, and compared to Monte Carlo, LHS, Halton, Sobol' and a rank-1 lattice. LAT is better than Monte Carlo and LHS as soon as the integrand has interactions. Sobol' remains more accurate at the powers of two, but its error is one to two orders of magnitude larger at other sample sizes. The accuracy of LAT does not depend on the sample size. Finally, the cells form a partition of the hypercube for any N. This allows one to refine the design locally, to search for an optimum by splitting cells and to sample non-rectangular regions.

We study false-alert control when screening for text generated by artificial intelligence (AI). The screening procedure selects document prefixes and detectors from observed evidence and may stop before exhausting its inspection budget. We give two finite-sample constructions under document-level exchangeability between human calibration documents and a new null document, with no restriction on dependence among tokens within a document. Construction A registers a finite family of prefix-detector scores and allocates a false-alert budget across their conformal ranks. A union bound protects any executed subset of that family. Construction B calibrates the complete-path maximum of a development-fixed adaptive policy. Each partial-path maximum is bounded by the complete maximum, so a terminal conformal rank protects early stopping without splitting the error budget. We prove marginal control of any false alert across the permitted inspection path and derive necessary calibration counts for rejection. We also state oracle testing, distribution-shift, and independent-audit bounds with their additional assumptions. Both constructions protect stopping within their specified scope; neither proof constructs an e-process or justifies multiplying conformal ranks. Detection power and computational savings remain questions for empirical evaluation.

The decision problem is whether valence among the users who write reviews of an app has changed enough to warrant investigation of a release, outage, or support issue. We propose a statistical specification for that decision, applied to written Google Play reviews collected under a declared locale, language, sort order, retrieval cap, and window rule. That endpoint returns a selected collection rather than a random sample, so the index describes review writers under the declared protocol and not all app users. What is new is not a sentiment classifier but an auditable interface: an ordinal threshold model is the preferred representation of star ratings, a calibrated two-signal generalized least-squares estimator is retained as an interpretable linear approximation when held-out labels support its conditional-unbiasedness assumption, the estimand is declared before the estimator, and every downstream module–uncertainty, diagnostics, and regularization–is tied to that declaration. The framework separates the latent average of a fixed collected set from the mean of a protocol-defined population of review writers, so measurement and between-review variation are not conflated, and treats helpfulness and recency weights as a different policy target. Shrinkage, cross-app empirical Bayes, histogram and dependence diagnostics, and dynamic smoothing are specified as separate modules with defaults and mandatory sensitivity analyses. A reference implementation, a declared 6000-review three-app collection, and a synthetic stress test exercise the displayed equations: fusion reduced error in clean and small-window scenarios, coordinated contamination caused severe bias and interval failure, and the collected sample reproduces the abstention rules. No real-review validation is claimed.

The standard dependence summaries used in biomarker studies – Pearson's r, Spearman's rho, Kendall's tau – take values in [-1, 1] with 0 indicating no linear or monotone association. Zero does not distinguish independence from non-monotone dependence, so the scale cannot represent threshold effects, heteroscedasticity, and tail shifts common in biomarker practice. We formulate independence testing as a composite Bernoulli likelihood ratio: at each threshold t, comparing the Bernoulli laws of 1(Y <= t) conditionally on X versus marginally, aggregated over cut points. The resulting coefficient xi_cut lies on [0, 1] with 0 iff X and Y are independent (under continuity of Y) and 1 iff Y is a measurable function of X. Fisher weighting arises at second order from the Bernoulli likelihood, and xi_cut equals twice the threshold-averaged mutual information between X and 1(Y <= t), giving a distribution-free lower bound on I(X; Y). A second-order expansion recovers the Fisher-weighted Dette-Siburg-Stoimenov measure, which coincides under continuity with Chatterjee's rank correlation. Estimation uses a Nadaraya-Watson plug-in with a max-over-grid bandwidth; inference is by exact permutation. In biomarker-motivated simulations T_cut substantially outperforms rank-based coefficients on W-shaped non-monotone and heteroscedastic alternatives. We illustrate on the Seattle cohort (n=70, ages 21-88) of the aging plasma proteome dataset, screening all 1,305 proteins for age dependence: under Benjamini-Hochberg control at q<0.05, T_cut rejects on 70 proteins, six of which are missed by Pearson, Spearman, and Chatterjee at the same FDR level.

Many studies expose units to several binary treatments at once and record several correlated outcomes, yet analysts usually estimate one treatment's average effect on one outcome at a time – discarding how treatments interact and how their effects vary across units. We introduce the factorial multivariate Bayesian causal forest, a Bayesian nonparametric model that decomposes the factorial response surface over the treatment lattice into a prognostic sum-of-trees plus one sum-of-trees per main or interaction effect, each with correlated multivariate leaf parameters, estimated jointly with coherent uncertainty. A single indicator-weighted kernel samples every component – the prognostic term is the special case whose indicator is unity – so the two-treatment model, the single-treatment multivariate causal forest, and a general order-r truncation are one and the same sampler. We give the ANOVA/Möbius identification of each estimand, an efficient Rcpp engine, and an interpretability layer using value-suppressing uncertainty maps. In simulations the average-effect estimators are unbiased with nominal coverage, consistent, and robust under misspecification, failing only under unmeasured confounding, which we flag. We illustrate the method on clinical, agricultural and economic data and on a deeper application to the NHANES survey. Software is provided as an R package.

Gaussian process prediction is a central tool in spatial statistics, but standard implementations require dense matrix operations that become prohibitive for large datasets. We propose a scale-adjusted compactly supported working correlation for Gaussian process prediction, using the generalized Wendland family with a support radius $\phi_n$ that is allowed to decrease with the sample size $n$. Under fixed-domain asymptotics with quasi-uniform designs, a fixed support radius does not yield asymptotic sparsity, whereas a shrinking support radius can make the covariance matrix sparse. We show that $\phi_n$ and the regularization parameter can be jointly chosen so that the resulting predictor preserves the optimal integrated mean squared prediction error rate while reducing the number of nonzero covariance entries and the associated sparse matrix-vector multiplication cost. We also establish an analogous rate-preserving sparsification result for kernel ridge regression. Simulations and an ERA5 temperature application show that the proposed method achieves competitive prediction accuracy while retaining the computational efficiency provided by sparse linear algebra.

Reliable high-dimensional variable selection requires scalable error-controlling methods. The Terminating-Random Experiments (T-Rex) selector estimates the false discovery rate (FDR) by aggregating early-terminated forward-selection paths in which predictors compete with synthetic dummies. We address two remaining challenges: i) predictor dependence can bias dummy-predictor competition; ii) computation is wasted on recomputing terms shared across experiments. Building on memory-efficient virtual dummies, which sequentially sample projections from their exact conditional law, we estimate the conditional Gaussian law of inactive predictors and map draws onto the remaining sphere radius. A shared lazy Gram cache computes the response product and each requested Gram column once across experiments. Simulations across three covariance structures show that uniform spherical dummies can exceed the target FDR, whereas the proposed method empirically controls FDR. Caching yields more than a sixfold speedup, and the method remains feasible with 100000 predictors, where competing FDR-controlling methods become computationally impractical.

Large language models (LLMs) as a judge enable scalable evaluation, but their judgments can be sensitive to response order and, even after removing such position effects, can still diverge systematically from human preferences.We introduce DIAL, a unified framework that combines abundant LLM comparisons with limited human comparisons to separate judge-specific position effects, learn shared structure in position-debiased LLM preferences, and adaptively calibrate that structure toward the human preference target. Theoretically, we study three aspects of DIAL: (i) identification of latent LLM preferences, position effects, and human calibration; (ii) adaptive estimation that balances LLM anchoring against limited human evidence; and (iii) fixed-weight uncertainty quantification for the calibrated human preference. Empirically, we evaluate position debiasing and human alignment separately in controlled simulations and on three human-preference benchmarks, showing that DIAL remains robust to unbalanced response order, achieves strong human-aligned rankings with limited labels, and adapts toward human evidence when LLM information is imperfect. Our real-data study collects over 410K judgments from 21 LLM judges in both display orders, providing a resource for future studies of LLM-judge bias, heterogeneity, and human alignment.

This paper develops prediction models for population-level survival curves observed over time and sampled from a heterogeneous mix of populations. We consider a discrete-time setting where each curve is only partially observed and forecasts of the remaining trajectory are needed for downstream decision making. Our approach recasts cross-population heterogeneity into a multivariate sampling model. We propose two forecasting models for partially observed curves: a full factor analysis model that extends a general factor representation to incorporate the partially observed survival curve, and a double PCA model. The methodology is motivated by demand forecasting in short-term rental markets, where market-level occupancy paths can be viewed as survival curves over the booking horizon and where forecasts of future occupancy feed into dynamic pricing algorithms. We apply the models to the newly released Wheelhouse dataset, which contains time series of market occupancy curves for 500 markets from 2017 to 2022. Model performance is assessed using the integrated quadratic distance, and we compare the proposed PCA-based methods to Holt's linear trend model across multiple forecast horizons. The results show that the proposed models yield accurate and stable forecasts of the remaining survival trajectory and generally outperform Holt's method, particularly at longer horizons.

LLM-as-a-judge has become the de facto standard for scalable, subjective evaluation, yet current leaderboards compensate for systematic measurement bias by running ever more comparisons, an approach that is both statistically unsound and computationally wasteful. The root cause is an incomplete measurement model, treating LLM judges as neutral, interchangeable instruments ignores documented biases like position bias, verbosity bias, judge severity, and self-enhancement, that no volume of additional data can eliminate. We propose a unified latent variable framework that jointly models pairwise and ordinal data while explicitly correcting for these confounders, recovering reliable rankings from substantially fewer comparisons. Because fitting this model costs negligible compute relative to a single round of LLM inference, bias correction is not only more statistically rigorous but also a more sustainable approach to trustworthy evaluation.

We study nonparametric multiple change-point detection for high-dimensional sequences, aiming to identify time points at which the underlying distribution changes. While many existing methods perform well when change-points are well-separated, their performance can deteriorate when structural breaks are densely clustered. To address this challenge, we propose gBottomup, a graph-based bottom-up framework for multiple change-point detection in high-dimensional settings. gBottomup constructs a hierarchical segmentation by proposing merges of adjacent segments and verifying them through an unmerge rule that combines absolute significance with relative local heterogeneity, thereby adaptively refining partitions and estimating the number of change-points. Simulation results demonstrate that gBottomup performs reliably across a range of structural configurations and is particularly effective in frequent change-point settings, where existing top-down procedures may lose sensitivity. Runtime experiments indicate favorable computational performance relative to graph-based top-down alternatives. We illustrate the proposed method through an analysis of a S&P 500 dataset.

Bipartite A/B tests are experiments in which treatment is randomised over one set of units, while outcomes are measured on another. For example, an online marketplace may test a new pricing algorithm on a random subset of items, while the outcome of interest, say, purchase satisfaction, is measured on customers, each of whom interacts with many items. Existing methods for analysing bipartite experiments assume that every randomisation unit is assigned to treatment or control. In practice, often only a subset participates: platforms cap rollout risk, reserve holdout groups, and split their population across concurrent tests. Ignoring the unassigned units biases estimation, while including them requires care. We construct an unbiased linear estimator for bipartite experiments with partial assignment. Observations must be reweighted not only by the (centred) share of treated connections among participating ones (the exposure) but also by the number of participating connections, each inverse-weighted by its participation propensity, so that units whose connections are well covered carry proportionally more weight. The resulting estimator, EARL (Exposure- and Allocation-Reweighted Linear), is unbiased, consistent, and asymptotically normal; we devise two asymptotic variance estimators and show that it has minimal variance in a natural class of linear estimators. EARL remains unbiased regardless of the experience the unassigned units receive, as long as they contribute to expected outcomes additively; in particular, they can be allocated to other, non-overlapping tests. Our theory is complemented with a simulation study on two public datasets, in which EARL attains up to six times lower error than the strongest existing baseline and, in some configurations, over an order of magnitude lower than Horvitz-Thompson-style alternatives.

Trimming methods for robust clusterwise regression discard a fixed fraction of the data. Too low a level breaks the fit; too generous a level can trim away a small group. We first propose a group-recovery step that can follow any trimming or flagging method: it searches the discarded units for a line, tests whether those near it form a peak rather than a band, and restores the line as a group when the likelihood of Gaussian groups plus uniform noise improves by a margin like that of the Bayesian information criterion. In simulations it repaired the failures of a generous TCLUST-REG level with unequal groups, but not with three or four groups. The second proposal, ESF (exact-subsample flagging), is a flagging procedure without a trimming level: it solves the clusterwise least-squares problem exactly on small subsamples, flags units far from the best fit, and draws later subsamples from the rest. Two constants stand in for the level: a subsample size, set from a lower bound on the smallest group proportion, and a cap on the flagged set. It is meant for data about whose contamination nothing is known: TCLUST-REG at the fixed level 0.30, followed by reweighting and the recovery step, was as accurate as ESF on average up to a fifth of outliers, and higher levels were more accurate beyond. The flagged fraction estimates the contamination only when the errors are close to Gaussian. On taxi fares with known tariffs, ESF found both in every sample.

Prediction models can perform poorly when the deployment population differs from the training population. Data from the target population would help, but individual-level target data may be inaccessible because of access restrictions or reporting conventions. We consider a multi-resolution setting in which individual-level data are available from a source population, while the target population is observed only through subgroup summaries. We propose SAGE, a one-step estimator that updates a source-trained predictor using a gradient estimated from these summaries. Motivated by diagnostics consistent with the random distribution shift model, we choose SAGE's step-size to account for both sampling and distributional uncertainty. Under this model, SAGE reduces mean asymptotic target excess risk relative to the source-trained predictor. We also show that, under the model, entropy-balancing weighted empirical risk minimization (EB), which reweights source observations to match the target summaries, is asymptotically equivalent to a full-step SAGE update. SAGE with the optimal step-size has asymptotic mean squared error no larger than that of EB. Across real-world datasets, SAGE generally improves on the source-trained predictor and one-step updates that ignore distribution shift, including when the random shift model only partially captures the observed shifts. Compared to EB, SAGE improves prediction more consistently across the sample size and shift settings studied.

Conformal prediction can lose coverage when the data distribution changes after deployment. We study adaptation using labeled source data and unlabeled target inputs, allowing both the input distribution and its relationship with outcomes to change. We use Exponential Tilt Reweighting Alignment (ExTRA), introduced for classification by Maity et al. (2023), to estimate structured distribution shifts. We compare using its estimated weights in conformal calibration with additionally tilting the source predictive distribution. Shared learned predictors, estimated weights, calibration samples, and test observations isolate the effect of tilting. Existing theory gives both procedures target coverage with true weights and a common coverage bound with estimated weights. Identification calculations and an analysis of how scoring interacts with weight estimation error help explain why their performance can nevertheless differ. In a synthetic regression setting where the assumed models match the data-generating process and target inputs are informative about the shift, tilting reduces mean set length by about $30\%$ relative to weighting alone, with both methods attaining coverage near nominal. Tilting can instead cause substantial coverage losses in synthetic classification and in regression when target inputs provide little information about the response shift. Real-data experiments also show no consistent benefit. Good coverage from weighted calibration alone does not ensure that adding predictive tilting will preserve coverage. Deciding when to apply this additional adjustment using only source labels and target inputs remains an open problem.

Multitask linear models can improve estimation and prediction by exploiting structure shared across responses, bridging taskwise fitting and complete pooling. In many applications, however, the objective is estimation in a data-limited target domain, while data-rich but heterogeneous source domains are available. Borrowing from these sources can improve efficiency but introduce bias. We develop a joint-sparse transfer-learning framework for high-dimensional multi-output regression that combines shared predictor structure across responses with source-target similarity. The framework yields two complementary estimators: a fused estimator that aggregates jointly fitted domain-specific coefficients and a target-based debiased estimator that adjusts for source-induced shifts. Our error bounds show how transfer increases the available information and sharing predictors across responses reduces selection costs. They also reveal a tradeoff: the fused estimator benefits from larger sources but may retain bias if source shifts point in similar directions, whereas debiasing trades this bias for additional estimation error governed by the smaller target sample. Comparison with a minimax lower bound identifies regimes where the bounds match up to logarithmic factors, where matching remains unresolved, and where projection onto a target-based convex set closes the gap. Simulations and an analysis of single-cell RNA and surface-protein profiles across cell types support the theory.

Neurons in the visual system are selective for specific spatial and temporal stimulus features, described by their receptive field. Estimating one means a coefficient per pixel per time bin from few trials – a high-dimensional problem requiring regularization. Sparse regularizers such as the LASSO handle the dimension but select pixels independently at each time point, with nothing to keep the region coherent in space or smooth in time; it can fragment or reorganize discontinuously even when the true response evolves smoothly, a failure since this evolving pattern is what a receptive-field estimate should capture. We formulate dynamic receptive-field estimation as a high-dimensional Bayesian problem: a Poisson model combining a Gaussian Markov random field in space with an autoregressive process in time, so the estimated field is smooth and coherent across space and time. On recordings from $155$ salamander retinal ganglion cells, fitting this model independently per neuron recovers a coherent surface, where a pixel-level Poisson-LASSO comparison instead returns a fragmented one. Summarizing each neuron's surface by its space-averaged temporal response and clustering these curves with a model-based functional-clustering procedure, BIC selects three balanced temporal-response phenotypes ($85$, $32$, $38$ neurons), against a degenerate grouping from clustering the raw surfaces. A simulation study with known ground truth confirms the same pattern, with the model beating an unregularized Poisson GLM, LASSO, and the elastic net on recovery and estimation accuracy, though LASSO controls false positives better. The per-neuron field identification, its contrast with LASSO, and the functional-clustering population typing constitute this paper's contribution.

In this work, we model the observed multi-view tensors by decomposing the underlying signal in each view into two components: (i) the shared component that captures common dynamics across all views, and (ii) the private component that accounts for view-wise unique variations. To decouple the shared and private components, we introduce a novel Tucker personalized subspace principal component analysis (TPS-PCA) approach for tensors, which admits a one-step closed-form solution and serves as an ideal surrogate for our extended tensor-version personalized PCA (TP-PCA), adapted from the seminal work by \cite{shi2024personalized}. The theoretical analysis reveals that the proposed TPS-PCA estimators reach the minimax lower bound in terms of view-wise tensor decoupling, whereas the TP-PCA estimators only achieve a rate of average decoupling error across views, which is still slower than that of the TPS-PCA estimators. Extensive numerical experiments are conducted on synthetic and real datasets, demonstrating the wide applicability of the proposed method in fields including power management, financial analysis, and activity recognition.

World models for epidemiology learn from records shaped by study designs, including assignment, sampling, measurement, and related processes. A model may therefore reconstruct observed trajectories while learning an intervention contrast that depends on how records were collected. We formalize study design as constraints on a world model latent-world interface, action mechanism, observation likelihood, and target readout. This distinguishes design-respecting models, which encode these constraints, from design-ignoring models, which fit selected records without relating assignment and observation to the intended intervention question. Using a large scale cluster-randomized test-negative trial, we hold the latent structure and fitting settings fixed and compare a design-ignoring case-count model with a design-respecting test-negative observation model. Both models achieve comparable factual reconstruction. Yet across 500 paired resampling experiments at different relative sampling intensity, median contrast changes are substantial for the design-ignoring model, versus merely marginal for the design-respecting model. The same qualitative separation holds when sampling also varies across clusters. Thus, factual reconstruction alone does not establish design alignment; testing whether fitted contrasts preserve design-implied observation-process invariances provides a sharper evaluation.

Spatial disease maps help public-health researchers identify geographic inequalities, but standard Bayesian smoothing can obscure localized disparities when neighboring communities have sharply different socioeconomic or behavioral profiles. Analysts therefore need to determine where smoothing should be interrupted and repeat that analysis as maps, adjacency structures, and outcomes change. We develop a covariate-informed Bayesian boundary model and an amortized posterior approximation trained across heterogeneous areal graphs. The model distinguishes local interruptions in smoothing from broader residual spatial dependence; the trained approximation handles maps with different numbers of regions. Simulations examine posterior calibration, boundary-probability recovery, replicated-data behavior, and MCSE-controlled agreement with prior-matched MCMC. In contrast to traditional approaches that analyze these data separately, we demonstrate the effectiveness of using a single trained deep learning network to analyze respiratory hospitalizations in Greater Glasgow; lung cancer incidence in California; and tracheal, bronchial, and lung cancer mortality in South Korea, comprising 58 to 241 regions. Selected boundary density is greatest in Glasgow and lowest in South Korea despite substantial residual spatial dependence in both, showing that local interruption and broader spatial persistence need not vary together. Across all three applications, edge-level boundary probabilities agree substantially with dataset-specific analyses, although posterior spread and thresholded boundary sets differ. These results support reusable Bayesian boundary analysis across the evaluated disease-map class and identify the validation needed before deployment to new applications.

We consider the problem of learning the underlying causal directed acyclic graph (DAG) structure corresponding to a structural equation model (SEM) with non-Gaussian errors. Motivated by an intentionally misspecified non-Gaussian SEM with all Laplace errors, we first introduce the mean absolute residual risk, defined over the space of all real matrices, and show that, asymptotically, the risk of the true weighted causal DAG matrix is strictly smaller than that of any other matrix. Nevertheless, to enhance generality and account for high-dimensional and finite-sample settings, we further incorporate row-specific sparsity penalties along with a soft DAG constraint to derive a continuous score function over the space of real matrices. Accordingly, we propose a score-based DAG learning method, named MARCEDES, formulated as an unconstrained score minimization problem, which can be efficiently solved using gradient-based optimization techniques, thereby circumventing the challenges associated with constrained optimization. Furthermore, we develop a computational algorithm to handle the non-smoothness of the score objective and to enable optimal tuning of row-specific sparsity penalties under a generalized Bayes framework. Finally, we demonstrate the efficiency and improved performance of the proposed method over existing approaches through an extensive simulation study.

Two-sample Mendelian randomization (MR) is a widely applied methodology in epidemiology. In two-sample MR, summary data (typically, regression coefficients and standard errors) quantifying the association between multiple genetic variants and the exposure and the outcome are used in an instrumental variable framework aimed at estimating the causal effect of the exposure on the outcome. Most two-sample MR methods were developed under data-generating models where the association of for each candidate genetic instrument with the exposure, as well as the causal effect of the exposure on the outcome, are constant in the additive scale. These assumptions are useful because they imply that, had all genetic variants been valid IVs, they would all estimate the same causal parameter - namely, the constant causal effect. We refer to this condition as summary-level homogeneity. However, these are rather strong homogeneity conditions which may raise concerns about the plausibility of these methods in practice. In this paper, we show that summary-level homogeneity is implied by the following conditions: the causal effect is additive linear, but not necessarily constant across, all strata of the population; and uncorrelatedness between heterogeneity in the causal effect and in the association between each genetic variant and the exposure. Under these conditions, typical two-sample MR methods can be interpreted as estimators of the average causal effect. These results clarify that point-identifying assumptions required for two-sample MR methods are weaker than previously anticipated, which contributes to their plausibility and interpretation in at least some practical applications.

Maximum likelihood estimation of Markov deterioration hazard models for infrastructure asset management is sensitive to the numerical conditioning of explanatory covariates, yet production pipelines often adopt a single scaling convention without systematic justification. We study four covariate scaling transformations—baseline max scaling, min-max scaling, z-score scaling, and Box-Cox transformation with a training-derived positivity shift—applied to an Exponential Hazard Markov (EHM) model estimated via L-BFGS-B on bridge inspection transition data.

We study the impact of conditional complier average causal effect (CCACE) estimation methods on the performance of subgroup discovery and heterogeneous causal effect estimation under imperfect compliance. Building on the Bayesian Causal Forest with Instrumental Variable (BCF-IV) (Bargagli-Stoffi et al. (2022)) method, we introduce a two-step, model-agnostic approach that allows any suitable machine learning method to be used for the CCACE estimation in the first step. Specifically, we implement two non-Bayesian tree-based methods, both using forest-based learners: DRRF-IV, a debiased transformed-outcome regression-forest approach, and a GRF-based IV adaptation of the generalized random forest framework (Athey et al., 2019). Through a simulation study, we assess the precision, bias, and the ability to correctly identify the underlying subgroup structure of the proposed methods relative to BCF-IV. The results show that the non-Bayesian methods perform competitively across the considered simulation settings, with performance improving for larger sample sizes and moderately large treatment effects, while reducing computational runtime. We apply our new methods by revisiting an empirical study that examines the effect of prompt admission to intensive care units (ICU) on 28-day mortality across 48 UK National Health Service hospitals. While previous work finds no significant overall treatment effect, we investigate whether subgroups of patients may benefit more from prompt ICU admission.

Martingale posteriors and related predictive resampling methods replace the likelihood–prior pair used within Bayesian inference with a predictive model for future observations. These methods are simple to implement and increasingly popular due to their computational efficiency, but little is known about their ability to accurately quantify uncertainty. In this work, we study the concentration and calibration properties of the martingale posterior for general functionals, and show that credible sets can systematically undercover if the predictive algorithms are not carefully tuned. We propose a simple remedy: the bagged martingale posterior. Rather than starting every predictive path from the observed sample, we start the paths from random bootstrap resamples of the data and then amalgamate the resulting draws. Critically, this approach incurs no additional {simulation} cost compared to standard predictive resampling, and delivers conservatively calibrated credible sets. This scalability enables application to a range of challenging examples, including sparse high-dimensional regression and nonparametric conditional quantile models.

Many real-world high-dimensional time series exhibit long-memory, but Gaussian graphical model testing in this regime remains understudied. We develop a direct, data-adaptive test statistic for assessing conditional independence in the graph structure of stationary Gaussian time series. We establish a finite-sample, Berry–Esseen type Gaussian approximation bound for the statistic, which applies to both short-memory and long-memory time series. The testing procedure is fully data-adaptive using block bootstrap method, on which we provide a finite-sample validity result including in the ultra-high-dimensional scenario, and can be extended to comparing graphical structures in two-sample tests. We also develop a consistency-empowered correction to the statistic and show that such tests attain asymptotic consistency in both size and power. Our proposed method is applied to a real-world fMRI data to understand functional connectivities within brain in different periods.

Benchy is a semantic language and execution engine for benchmarking AI programs. A benchmark is completely specified by a program, a scoring function, and a dataset, B=(P,S,D), and is separate from the AI-system taking it; a run binds the two, R=(B,AI). Benchmarks are authored as canonical YAML in which each semantic concept has one valid syntax, classified by a shared task/domain/language ontology, and deterministically compiled into a canonical JSON intermediate representation that the engine executes. Compilation changes representation, not meaning: it does not repair invalid definitions or inject hidden defaults. Programs use fixed schemas of named input and output fields, the leaf output fields are the scoring dimensions, and the engine exposes one universal runtime contract — a named-field input object in, a named-field output object out — to which external AI-systems adapt at the boundary, so integration mechanics never propagate into benchmark semantics. This paper gives the semantic object model, the ontology and task-to-program validation rule, the scoring and failure semantics, the compilation and execution architecture, and the scope of the current language. An appendix fixes the normative engineering contract for the first engine implementation.

Online experiments must often be evaluated before long-term outcomes mature. Under rolling enrollment, these outcomes are observed only for early enrollees, while short-term surrogates are available for everyone. We compare seven estimators across eleven data-generating processes, spanning partial mediation, drift, outcome sparsity, and enrollment-time labeling, with up to $R = 2,000$ replications over more than 500 method-by-scenario cells. We find a sharp robustness-efficiency tradeoff: the surrogate index delivers large efficiency gains when surrogacy holds but its coverage collapses under violations, while PPI-family methods stay asymptotically valid under random labeling at smaller gains. We give the finite-sample variance of PPI++ in the all-units parameterization for a fixed predictor, a joint asymptotic distribution for the two estimators under cross-fitting, and a Hausman-type estimator-disagreement diagnostic, then quantify the detection-damage gap: in the partial-mediation design, where we locate both edges, a band of violations destroys surrogate-index coverage yet is too small to detect on most datasets. On the 64,000-customer Hillstrom experiment the diagnostic rarely flags a violation that biases the surrogate index, and a Cauchy-kernel hybrid of the two estimators inherits 10.7% relative bias; on the 14-million-user Criteo experiment the violation is detected, and subsampling traces detection turning on with scale as damage persists. PPI++ has limits: at a rare-conversion $n = 30,000$ Criteo subsample its empirical coverage is 85.0%. We recommend prespecifying PPI++ with the exact variance as the primary analysis under random labeling and adequate labeled outcome counts, reading the diagnostic as a warning, not a certificate, and treating the surrogate index as a sensitivity analysis.

Breast cancer risk assessment and prediction models are very important tools for the clinical management of healthy women carrying pathogenic variants in known cancer genes. We propose in this paper a comprehensive prediction model for the personalized clinical management of women with pathogenic variants in BRCA1. This model is able to estimate the risk of BC accounting for the full history of this cancer within the family, the personal and family history of risk-reducing salpingo-oophorectomy as well as the exact age of this intervention among family members. The effect of oophorectomy on the risk of developing breast cancer is evaluated by estimating and conducting inference about the regression coefficients in a Cox model with time-varying covariates. We model the within-family dependence in the ages at onset of cancer using a Gaussian copula whose correlation matrix accommodates the different pairwise family relationships. We develop an iterative algorithm to estimate the parameters of the considered model in the presence of a selection bias. We evaluated the performance of the proposed cancer risk estimation method from family data by simulations and illustrated its use through an application to the breast cancer family registry.

This paper presents an end-to-end causal machine learning (ML) pipeline designed for real-world applications with continuous treatments. The proposed framework consists of six sequential steps: dimensionality reduction, causal identification, positivity assumption violation handling, estimation, refutation and evaluation, and policy optimization. We introduce practical contributions not currently available in existing causal ML toolkits, specifically: (1) a method for detecting and quantifying positivity violations in continuous treatment settings (2) a novel, scalable two-stage dimensionality reduction framework tailored for causal inference with high-dimensional data; (3) the adaptation of sensitivity analysis and estimation methods originally designed for binary treatments to the continuous treatment space and (4) an end-to-end integration of these components into a modular, reproducible workflow. These innovations address real-world challenges in causal inference that are often not covered in theoretical frameworks but frequently encountered in industrial applications. The methodology is validated with a synthetic dataset inspired in a real-world financial debt collection use case, however its design can be applied to analogous problems across different industries. Results demonstrate that the proposed methodology offers a more computationally efficient approach and produces less biased estimates compared to standard methods for problems with continuous treatment and high-dimensional data. A fully functional GitHub repository with documented code and numbered notebooks is made available ensuring reproducibility and practical implementation. The pipeline presented is intended to contribute to closing the gap between academic approaches and practical application in industry contexts where causal ML can be highly beneficial such as the financial sector.