Fits isotonic distributional regression (IDR) to a training dataset.
idr(y, X, y_observed = NULL, weights = NULL, decreasing = FALSE,
groups = setNames(rep(1, ncol(X)), colnames(X)), orders = c("comp" = 1),
stoch = "sd", pars = list(verbose = FALSE, eps_abs = 1e-5,
eps_rel = 1e-5, max_iter = 10000L), progress = TRUE)An object of class "idrfit" containing the following
components:
Xthe training covariates as provided, one row per observation (in input order, including duplicated rows).
ynumeric vector of the training responses.
cdfmatrix containing the estimated CDFs, one CDF per row,
evaluated at response_unique (see next point). The CDF in the
ith row corresponds to the estimated conditional distribution of the
response given the covariates values in X[i,].
weightsthe observation weights as provided (NULL if
none were given).
response_uniquethe thresholds at which the CDFs in
cdf are evaluated. The entries in cdf[,j] are the conditional
CDFs evaluated at response_unique[j].
groups, ordersthe groups and orders used for estimation.
diagnosticdiagnostics of the CDF estimation. For univariate
fits (total order) this is list(epsilon = ), a bound on the
precision of the CDF estimation (the maximal downwards-step in the CDF
that has been detected). For multivariate fits (partial order) this is
list(precision = , convergence_fraction = ), where
convergence_fraction is the fraction of CDF estimations that
converged before hitting the iteration limit max_iter. Decrease the
parameters eps_abs and/or eps_rel or increase
max_iter in pars to improve the precision.
numeric vector (the response variable).
data frame of numeric or ordered factor variables (the regression covariates).
vector of indicators (TRUE or 1 for observed, FALSE or 0
for right-censored). At least one observation must be uncensored. Default
is all observed (rep(TRUE, length(y))).
vector of finite, non-negative weights (same length as y), at least one of which must be positive; observations with zero weight are dropped from the fit. Default is all weights equal to one. Weights are processed in single precision; it is up to the caller to avoid extreme imbalance (as a rule of thumb, no weight below ~1e-7 of the total weight).
boolean indicating whether y decreases with X
(by default, it increases with X).
named vector of length ncol(X) denoting groups of
variables that are to be ordered with the same order (see 'Details'). Only
relevant if X contains more than one variable. The same names as in
X should be used.
named vector giving for each group in groups the order
that will be applied to this group. Only relevant if X contains more
than one variable. The names of orders give the order, the entries
give the group labels. Available options: "comp" for componentwise
order, "sd" for stochastic dominance, "icx" for increasing
convex order (see 'Details). Default is "comp" for all variables.
The "sd" and "icx" orders can only be used with numeric
variables, but not with ordered factors.
stochastic order constraint used for estimation. Default is
"sd" for first order stochastic dominance. Use "hazard" for
hazard rate order (experimental).
parameters for quadratic programming optimization (only relevant
if X has more than one column), a list with options "verbose" T / F
(verbosity of solver), "eps_abs" positive float, "eps_rel" positive float,
"max_iter" positive integer.
display a progress bar while fitting (TRUE,
FALSE or 1, 0). Default is TRUE; the bar is
written to stderr and is best viewed in an interactive R session.
This function computes the isotonic distributional regression (IDR)
of a response y on on one or more covariates X. IDR estimates
the cumulative distribution function (CDF) of y conditional on
X by monotone regression, assuming that y is more likely to
take higher values, as X increases. Formally, IDR assumes that the
conditional CDF \(F_{y | X = x}(z)\) at each fixed threshold z
decreases, as x increases, or equivalently, that the exceedance
probabilities for any threshold z \(P(y > z | X = x)\) increase
with x.
The conditional CDFs are estimated at each threshold in unique(y).
This is the set where the CDFs may have jumps. If X contains more
than one variable, the CDFs are estimated by solving
length(unique(y)) quadratic programs with a built-in operator
splitting solver (see references). This might take a while if the training
dataset is large.
Use the argument groups to group exchangeable covariates.
Exchangeable covariates are indistinguishable except from the order in
which they are labelled (e.g. ensemble weather forecasts, repeated
measurements under the same measurement conditions).
The following orders are available to perform the monotone regression in IDR:
Componentwise order ("comp"): A covariate
vector x1 is greater than x2 if x1[i] >= x2[i] holds
for all components i. This is the standard order used in
multivariate monotone regression and should not be used for
exchangeable variables (e.g. perturbed ensemble forecasts).
Stochastic dominance ("sd"): x1 is greater than x2 in
the stochastic order, if the (empirical) distribution of the elements of
x1 is greater than the distribution of the elements of x2 (in
first order) stochastic dominance. The "sd" order is invariant under
permutations of the grouped variables and therefore suitable for
exchangeable covariables.
Increasing convex order ("icx"):
The "icx" order can be used for groups of exchangeable variables. It
should be used if the variables have increasing variability, when their
mean increases (e.g. precipitation forecasts or other variables with
right-skewed distributions). More precisely, "icx" uses the
increasing convex stochastic order on the empirical distributions of the
grouped variables.
Henzi, A., Moesching, A. & Duembgen, L. Accelerating the Pool-Adjacent-Violators Algorithm for Isotonic Distributional Regression. Methodol Comput Appl Probab (2022). https://doi.org/10.1007/s11009-022-09937-2
Bladt, M., Henzi, A., van den Heuvel, B. and Ziegel, J. (2026). Survival Isotonic Distributional Regression. arXiv:2608.02914. https://doi.org/10.48550/arXiv.2608.02914
The S3 method predict.idrfit for predictions based on
an IDR fit.
data("rain")
## Fit IDR to data of 185 days using componentwise order on HRES and CTR and
## increasing convex order on perturbed ensemble forecasts (P1, P2, ..., P50)
varNames <- c("HRES", "CTR", paste0("P", 1:50))
X <- rain[1:185, varNames]
y <- rain[1:185, "obs"]
## HRES and CTR are group '1', with componentwise order "comp", perturbed
## forecasts P1, ..., P50 are group '2', with "icx" order
groups <- setNames(c(1, 1, rep(2, 50)), varNames)
orders <- c("comp" = 1, "icx" = 2)
fit <- idr(y = y, X = X, orders = orders, groups = groups)
fit
Run the code above in your browser using DataLab