Package {geosmooth}

geosmooth package

Description

Geometric smoothing and conditional expectation methods for coordinate data, point clouds, and weighted graphs. The package provides local polynomial, model-averaged, trend-filtering, graph low-pass, occupation-density, and Hessian-energy methods.

Start here

Use browseVignettes("geosmooth") to open the installed HTML guides, including vignette("function-guide", package = "geosmooth"). The online guide links to the full reference. The task links below work in installed help without an internet connection.

Finding a function

Read vignette("function-guide", package = "geosmooth") for a task-based guide to all exported functions, prediction limits, and examples.

Local regression

Start with fit.lps or fit.malps.

Lifting trend filters

Use fit.lpl.tf or fit.slpl.tf for residual penalties with optional prediction synchronization.

Univariate trend filters

Use fit.gaussian.trend.filtering for continuous responses or fit.logistic.trend.filtering for binary responses. Assess association and prediction with assess.gaussian.trend.filtering or assess.logistic.trend.filtering, respectively.

Graph smoothing

Use fit.metric.graph.lowpass for spectral smoothing or fit.graph.trend.filtering for graph-difference penalties.

Hessian penalties

Use fit.ssrhe.hessian.regression for quadratic energy or fit.ssrhe.hessian.l1.regression for an absolute-value penalty.

Probability mass

Use fit.density for support weights, fit.subject.od for repeated visits, and normalize.density to convert an existing fitted field.

Synthetic examples

Combine components with synthetic.spec and draw data with materialize.synthetic.

Author(s)

Maintainer: Pawel Gajer pgajer@gmail.com

Authors:

• Pawel Gajer pgajer@gmail.com

Other contributors:

• Gael Guennebaud (Eigen) [contributor]

• Benoit Jacob (Eigen) [contributor]

• Authors of Eigen (Authorship and copyright in the included Eigen library) [copyright holder]

• Yixuan Qiu (Spectra) [contributor, copyright holder]

• Contributors to Spectra (Copyright in the included Spectra library) [copyright holder]

• Sunil Arya (ANN) [contributor, copyright holder]

• David M. Mount (ANN) [contributor, copyright holder]

• University of Maryland (ANN) [copyright holder]

See Also

vignette("function-guide", package = "geosmooth"), vignette("synthetic-datasets", package = "geosmooth")

Examples

X <- matrix(seq(0, 1, length.out = 30), ncol = 1)
y <- sin(2 * pi * X[, 1])
fit <- fit.lps(X, y, support.grid = 8L, degree.grid = 1L,
               kernel.grid = "tricube", foldid = rep(1:3, 10))
head(predict(fit))

Apply a Metric Graph Low-Pass Filter Path

Description

Applies every requested low-pass parameter to one response or a matrix of responses using a reusable spectral basis.

Usage

apply.metric.graph.lowpass.path(
  basis,
  y,
  eta.grid,
  filter.type = c("heat_kernel", "tikhonov", "cubic_spline", "gaussian", "exponential",
    "butterworth"),
  block.size = NULL,
  truncation.tol = 1e-04,
  unresolved.action = c("warn", "error", "allow"),
  exact.zero = TRUE,
  numerical.policy = NULL
)

Arguments

Value

A list of class "metric.graph.lowpass.path". For one response, fitted.values is an N by J matrix. For multiple responses it is an N by J by S array.

Examples

adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
lengths <- list(1, c(1, 1), c(1, 1), 1)
basis <- metric.graph.lowpass.basis(adj, lengths, n.eigenpairs = 4L,
                                    eigen.solver = "dense")
apply.metric.graph.lowpass.path(basis, 1:4, eta.grid = c(0, 0.1, 1))

Summarize a synthetic dataset as one metadata row

Description

Summarize a synthetic dataset as one metadata row

Usage

## S3 method for class 'synthetic_dataset'
as.data.frame(x, row.names = NULL, optional = FALSE, ...)

Arguments

Value

A one-row data frame.

Examples

x <- materialize.synthetic(synthetic.registry.spec("G1"), n = 20L, seed = 1L)
as.data.frame(x)

Test Constancy and Assess Gaussian Trend-Filtering Prediction

Description

Assess a continuous outcome using order-2 trend filtering, retuned outcome permutations and a separate nested-CV comparison with a constant mean.

Usage

assess.gaussian.trend.filtering(
  object,
  permutations = 999L,
  inner.folds = 5L,
  outer.folds = 5L,
  seed = 1L,
  predictive = TRUE,
  warn = TRUE,
  additional.tests = character()
)

## S3 method for class 'gaussian.trend.filtering.assessment'
print(x, ...)

Arguments

Details

The randomization null is an outcome distribution that is unchanged across the fixed predictor locations. Independent outcomes must be exchangeable under the null; a constant mean plus identically distributed independent errors is a sufficient model. Normality is unnecessary. This is not a test of constant mean under arbitrary heteroskedasticity. Repeated participants, serial dependence, covariate adjustment and unequal weights require other null procedures and are rejected where detectable. Independence and exchangeability cannot be checked from a fit.

Fresh balanced folds are generated independently of outcomes; stored fit selections and fold labels are ignored. The explicit lambda path is tuned by CV on the full data and again on every shuffled outcome vector, always using the same tuning folds. The primary statistic is max(0, sum((y-mean(y))^2) - sum((y-f.hat(x))^2)), the nonnegative reduction in residual sum of squares, in squared response units. It is not a penalized-objective difference or a chi-squared statistic. Constant outcomes have zero statistics and permutation p-values of one.

With B requested permutations the inclusive plus-one p-value is (1 + sum(permuted.statistic >= observed.statistic))/(B + 1), allowing a conservative numerical tie tolerance of 1e-10 * max(1, abs(observed.statistic)). Retuning is part of each randomization. Failed candidates are ineligible, and a replicate fails if none is eligible. Failed replicates are retained without retries or denominator changes. Any failure suppresses the point p-value; lower and upper bounds count missing statistics as nonexceedances or exceedances. Observed-fit failure leaves both bounds unavailable. These are failure bounds, not confidence intervals.

Absolute deviation is mean(abs(f.hat(x)-mean(y))), averaging over observations with repeated predictors counted at their observed frequency. Total variation sums absolute changes along the continuous fitted curve over ⁠[min(x),max(x)]⁠, including every knot and interior quadratic turning point. Both are in response units. Total variation is not the third- difference penalty. These describe the newly tuned observed fit; they are not bias-corrected population effects or uncertainty intervals. Optional tests reuse the same fitted permutations, with per-statistic failure bounds and Holm adjustment across squared error and the requested tests. An incomplete family has no Holm p-values. This adjustment does not cover additional predictors or outcomes. Choose tests and penalties in advance.

Nested CV inherits the input fit's cv.penalty.scale. Under "training_mass", outer training fits scale the reference lambda grid and ridge penalties by their sample-size fraction; inner folds apply their own training fraction. Full-data observed and permutation refits retain the reference grid. Older saved objects without this field use "none" to preserve their earlier convention. Nested CV tunes using only each outer training set. Its constant baseline is that set's sample mean. Pooled mean squared error (mse) and mean absolute error (mae) use the same held-out rows; improvements are constant minus smoother loss. Positive values favor trend filtering. The baseline remains the training mean for both scores. Boundary predictions use the fitter's constant-boundary policy and are flagged. Incomplete predictions leave overall smoother scores unavailable. No standard errors are computed from overlapping CV folds. Prediction scores assess training at the outer-training sample size, and need not improve when the global association test rejects.

Runs are sequential, with the same fitting cost and failure policy as assess.logistic.trend.filtering. The test does not establish nonlinearity, causality or a useful predictive effect.

Value

A gaussian.trend.filtering.assessment with the same evidence structure as assess.logistic.trend.filtering, replacing null.probability with null.mean and the primary likelihood test with squared.error. Scalar statistic, p.value, p.value.bounds, status and the permutations table refer to squared-error improvement. fit is the newly tuned observed fit. curve.summaries, tests, and permutation.statistics report curve summaries and requested tests. Fold assignments, seeds, candidate diagnostics and failure messages are retained. Optional prediction includes row-level mean and constant.mean predictions, losses, fold diagnostics, model scores and baseline-minus-smoother improvements. No confidence intervals are supplied.

References

Phipson and Smyth (2010). Permutation p-values should never be zero. doi:10.2202/1544-6115.1585

See Also

fit.gaussian.trend.filtering, assess.logistic.trend.filtering

Examples

x <- seq(0, 1, length.out = 40)
y <- sin(2*pi*x) + 0.2*cos(seq_along(x))
fit <- fit.gaussian.trend.filtering(x, y, c(Inf, 0.01), backend = "banded")
# Nine permutations illustrate the API; use more for inference.
assessment <- assess.gaussian.trend.filtering(fit, permutations = 9,
    inner.folds = 3, outer.folds = 3, seed = 42,
    additional.tests = c("absolute.deviation", "total.variation"))
assessment$tests
assessment$prediction$improvement

Test Constancy and Assess Logistic Trend-Filtering Prediction

Description

Assess whether an order-2 logistic trend filter detects a nonconstant binary probability, and separately estimate predictive improvement over a constant. This is an additional assessment, not an option that changes the fitted objective or supplies confidence bands.

Usage

assess.logistic.trend.filtering(
  object,
  permutations = 999L,
  inner.folds = 5L,
  outer.folds = 5L,
  seed = 1L,
  predictive = TRUE,
  warn = TRUE,
  additional.tests = character()
)

## S3 method for class 'logistic.trend.filtering.assessment'
print(x, ...)

Arguments

Details

The null is a common unknown Bernoulli probability at all observed predictor locations. The constant estimate is the sample event rate, not a uniformly weighted integral over the continuous predictor domain. Independent observations must be exchangeable under this null. Repeated participants, serial dependence, confounding adjustment and grouped binomial counts require other null procedures and are not supported here. Independence cannot be established from a fitted object.

The input object's selected penalty and fold labels are not reused: those labels may depend on outcomes. A new outcome-independent fold assignment selects a penalty from the supplied path, including a full-data refit. Every permutation repeats this complete selection procedure on the same folds with the same solver settings. The statistic is the nonnegative part of twice the ordinary binomial log-likelihood improvement over the constant. It is not a penalized-objective difference and is not compared with a chi-squared distribution. Lambda order determines selection ties. Changing the path after inspecting the assessment invalidates the stated test unless that additional selection is included in the randomization.

Absolute deviation is mean(abs(p.hat(x) - mean(y))), evaluated at all observed predictor locations, including repeats. It averages over sampled observations, not uniformly over the predictor domain. A value of 0.04 means an average absolute departure of four percentage points. Total variation is the sum of absolute probability changes along the fitted curve over ⁠[min(x), max(x)]⁠. It includes every knot and interior turning point of the piecewise quadratic logit, up to floating-point precision; it is not a coarse-grid approximation or the trend-filtering penalty. For a monotone curve it equals the endpoint probability difference in absolute value; an oscillating curve can have total variation above one. These are descriptive summaries of the newly tuned fit, not uncertainty intervals or bias-corrected population effect estimates.

Optional tests use the same selected fits, folds and permutations as the likelihood test, with large values indicating departure from constancy. No additional model fits or random draws are needed. Each statistic has its own inclusive comparison tolerance and failure bounds. tests reports individual p-values and Holm-adjusted p-values across likelihood and the requested additional tests. Holm adjustment remains valid under dependence; it does not adjust across other predictors or outcomes. If any requested test is incomplete, all Holm values are unavailable. Choosing the smallest unadjusted p-value after inspecting several tests does not preserve the stated test level. The legacy scalar fields always describe the unadjusted likelihood test, regardless of additional.tests.

Outcomes are shuffled uniformly with the predictor fixed, preserving the total event count. With B successful permutations the Monte Carlo p-value is (1 + sum(permuted.statistic >= observed.statistic))/(B + 1). Under exchangeability this is a valid conditional randomization test, subject to Monte Carlo discreteness. For numerical ties, comparisons allow a conservative tolerance of 1e-10 * max(1, abs(observed.statistic)). A failed candidate is ineligible under the fitter's usual CV rule. If no candidate is eligible, that replicate fails. Failed replicates are never dropped or retried: the point p-value is NA and p.value.bounds gives the smallest/largest possible Monte Carlo p-value counting each missing statistic as a nonexceedance/exceedance. These bounds are not a confidence interval. An observed-fit failure leaves both bounds unavailable.

Predictive assessment uses separate outer folds. Each outer training set inherits the input fit's cv.penalty.scale. Under "training_mass", its lambda grid and ridge penalties are multiplied by the outer/full sample-size ratio, and its inner folds apply their own training ratio. Observed and permuted full-data refits use the same unscaled reference grid. Older saved fits without this field retain "none" for reproduction. Each outer training set tunes the path using only its inner folds and predicts its outer test set. The constant predictor is that outer training set's event rate. Both use identical test observations. Pooled mean natural-log loss and Brier score are reported; improvements are constant loss minus trend-filtering loss, so positive values favor the smoother. These assess training on the outer training-set size, not the full sample. No standard error is computed from overlapping CV folds. Incomplete smoother predictions leave its overall scores and improvements unavailable rather than scoring a selected subset.

Computation is sequential. A complete run uses approximately (permutations + 1) * (inner.folds + 1) native path fits for the test, plus outer.folds * (inner.folds + 1) for prediction. Per-replicate diagnostics are retained, but full permuted fitted objects are not. Rejecting constancy is evidence of association, not necessarily nonlinearity, predictive usefulness or a causal effect.

Value

A logistic.trend.filtering.assessment list containing statistic, p.value, p.value.bounds, status, null.probability, the newly tuned observed fit, permutation statistics and diagnostic records, optional nested-CV prediction, realized folds, seeds and solver settings. curve.summaries is a named vector with absolute.deviation and total.variation; curve.summary.details records their conventions and any calculation failure. tests has one row per requested test, including likelihood, with statistics, status, completed/failed permutation counts, exceedances, comparison tolerances, raw p-values, lower/upper failure bounds and p.value.holm. permutation.statistics contains a replicate column and one column per tested statistic. The original permutations table retains its likelihood-only schema. prediction$scores contains model losses; prediction$improvement contains baseline-minus-smoother gains; prediction$predictions and prediction$diagnostics retain row-level and fold-level evidence.

References

Ojala and Garriga (2010). Permutation Tests for Studying Classifier Performance. https://www.jmlr.org/papers/v11/ojala10a.html

Phipson and Smyth (2010). Permutation p-values should never be zero. doi:10.2202/1544-6115.1585

See Also

fit.logistic.trend.filtering

Examples

x <- seq(0, 1, length.out = 60)
y <- rep(c(0, 1, 0, 1, 1), 12)
fit <- fit.logistic.trend.filtering(x, y, lambda = c(Inf, 0.01),
    backend = "banded", curve.ridge = 0.1)
# Only nine permutations for a quick illustration; use more for inference.
assessment <- assess.logistic.trend.filtering(fit, permutations = 9,
    inner.folds = 3, outer.folds = 3, seed = 42)
print(assessment)
assessment$prediction$improvement
assessment$curve.summaries
# Optional tests reuse the fits within an assessment run:
shape <- assess.logistic.trend.filtering(fit, permutations = 9,
    inner.folds = 3, seed = 42, predictive = FALSE,
    additional.tests = c("absolute.deviation", "total.variation"))
shape$tests

Bootstrap Uncertainty For A Fixed-Selection MALPS Fit

Description

Reuses the supports, local charts, prediction supports, averaging weights, selected degree, and local-solver controls from an existing fit.malps object, then repeatedly calls refit with bootstrap-style case weights. This is a fixed-selection uncertainty diagnostic: it measures refit variability conditional on the fitted MALPS support profile and does not rerun cross-validation or rebuild supports in each replicate.

Usage

bootstrap.malps(
  object,
  B = 200L,
  weight.type = c("bayesian", "multinomial"),
  y = NULL,
  probs = NULL,
  conf.level = 0.95,
  seed = NULL,
  max.failures = max(10L, B),
  keep.weights = FALSE,
  verbose = FALSE,
  ...
)

Arguments

Details

Two weight generators are available. "bayesian" draws positive exponential weights and rescales them to have mean one. This is the default because it preserves positive local support weights and therefore avoids many local-support degeneracies. "multinomial" draws ordinary nonparametric bootstrap counts with total count nrow(object$X); this can fail when a replicate assigns zero weight to all observations in a local support, and such failures are recorded.

Value

A list of class "malps_bootstrap" containing:

• replicate.fitted.values: an n \times B matrix of successful replicate fitted values;

• replicate.weights: an optional n \times B matrix of bootstrap weights when keep.weights = TRUE;

• failures: a data frame with columns attempt, failure, and message;

• summary: pointwise fitted values, bootstrap summaries, and empirical lower/upper intervals;

• bootstrap metadata such as B.completed, n.failures, conf.level, and weight.type.

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
fit <- fit.malps(X, X[, 1]^2, degree = 1L,
                 support.type = "knn", support.size = 8L)
bootstrap.malps(fit, B = 5L, seed = 1L)

Compare Masked Graph Fields with Independent Graph Workers

Description

Runs one candidate graph per process, with dynamically scheduled jobs and no nested fold parallelism. On Unix, forked workers inherit the loaded package sources. Run from Rscript, not a GUI or multithreaded parent process.

Usage

compare.masked.graph.fields(operators, studies, workers = 1L, ...)

Arguments

Value

Named graph results, each containing named study results or an explicit error record. Retains per-observation predictions for common-support scoring.

Examples

op <- metric.graph.lowpass.operator(list(2L, c(1L, 3L), 2L),
                                   list(1, c(1, 1), 1))
studies <- list(example = list(vertex = rep(c(1, 3), 6), y = 1:12,
                               group = letters[1:12]))
comparison <- compare.masked.graph.fields(list(path = op), studies,
    workers = 1, lambdas = c(0.01, 1), outer.folds = 3, inner.folds = 2)
comparison$path$example$coverage

Compare two synthetic datasets scientifically

Description

Compare two synthetic datasets scientifically

Usage

compare.synthetic.dataset(x, y, tolerance = c(1e-12, 1e-10))

Arguments

Value

A list with equality and mismatch details.

Examples

spec <- synthetic.registry.spec("G1")
x <- materialize.synthetic(spec, n = 20L, seed = 1L)
y <- materialize.synthetic(spec, n = 20L, seed = 1L)
compare.synthetic.dataset(x, y)

Participant-Grouped Nested Cross-Validation for a Masked Graph Field

Description

Select smoothing within outer training folds, predict held-out groups, and finally tune and refit on all observations. Graph construction is label-free and fixed (transductive prediction). Unsupported predictions remain NA.

Usage

cv.masked.graph.field(
  operator,
  vertex,
  y,
  group,
  lambdas = 10^seq(-6, 2, length.out = 17),
  outer.folds = 5L,
  inner.folds = 4L,
  seed = 1L,
  fold.id = NULL,
  progress = NULL
)

Arguments

Value

Outer predictions, baseline predictions, scores on supported cases, coverage, tuning records, fold assignments and final full-data field. Scores do not account for selecting a graph after comparing graphs.

Examples

op <- metric.graph.lowpass.operator(list(2L, c(1L, 3L), 2L),
                                   list(1, c(1, 1), 1))
cv <- cv.masked.graph.field(op, rep(c(1, 3), 6), 1:12, letters[1:12],
    lambdas = c(0.01, 1), outer.folds = 3, inner.folds = 2)
cv$coverage
cv$final$fitted

Fit a Chart-Kernel Smoother

Description

Fits a local chart-kernel field by evaluating a Nadaraya–Watson style smoother at each evaluation point. The function is a general fitted-field model, not an occupation-density wrapper: density workflows should call normalize.density on the returned fit.

Usage

fit.chart.kernel(
  X,
  y,
  X.eval = NULL,
  support.size = min(15L, nrow(X)),
  kernel = c("gaussian", "tricube", "epanechnikov", "triangular"),
  bandwidth.multiplier = 1,
  support.grid = NULL,
  kernel.grid = NULL,
  bandwidth.multiplier.grid = NULL,
  foldid = NULL,
  cv.folds = 5L,
  cv.seed = 1L,
  coordinate.method = c("coordinates", "local.pca"),
  chart.dim = NULL,
  chart.dim.grid = NULL,
  selection.strategy = c("grid", "sparse_kd", "plateau_kd"),
  chart.dim.max = NULL,
  geometry.margin = 0L,
  auto.chart.support.metric = c("coordinates", "operator", "both"),
  auto.chart.selection.metric = c("coordinates", "operator"),
  quadrature.weights = NULL,
  denominator.floor = sqrt(.Machine$double.eps),
  return.details = TRUE
)

Arguments

Details

For an evaluation point x_u, the method chooses a local support U_u, builds either centered ambient coordinates or a local PCA chart, and computes

\widehat f(x_u) = \frac{\sum_{i\in U_u} y_i K_h(z_{ui})} {\sum_{i\in U_u} q_i K_h(z_{ui})}.

Here z_{ui} is the local chart coordinate of x_i-x_u, q_i is an optional quadrature weight, and K_h is the selected kernel.

Value

A list with class "chart_kernel" containing fitted.values, source/evaluation supports, selected controls, and denominator diagnostics.

Examples

X <- matrix(seq(0, 1, length.out = 21), ncol = 1)
y <- exp(-20 * (X[, 1] - 0.5)^2)
fit.chart.kernel(X, y, support.size = 7L, kernel = "gaussian")

Fit a Density Estimator

Description

Dispatches to a dedicated density estimator. Subject-occupation density estimation is one application: construct a sparse mass vector over a fixed support set and call this generic density layer.

Usage

fit.density(
  X,
  weights = NULL,
  method = c("empirical", "graph_random_walk"),
  graph = NULL,
  graph.control = list(),
  density.control = list(),
  return.details = TRUE,
  ...
)

Arguments

Value

A list of class "density_fit" with fields method.id, status, rho, empirical.rho, fitted.raw, theta, accounting, smoothness, timing, diagnostics, and warnings.

Examples

X <- matrix(seq(0, 1, length.out = 6), ncol = 1)
fit.density(X, weights = c(1, 0, 2, 0, 0, 1), method = "empirical")

Fit Order-2 Gaussian Trend Filtering

Description

Fit a continuous response with a piecewise-quadratic curve and an L1 penalty on third divided differences. Dense and banded native C++ solvers share the order-2 geometry of fit.logistic.trend.filtering.

Usage

fit.gaussian.trend.filtering(
  x,
  y,
  lambda,
  weights = NULL,
  covariates = NULL,
  offset = NULL,
  nuisance.ridge = 0,
  curve.ridge = 0,
  lambda.selection = c("fixed", "cv"),
  foldid = NULL,
  control = list(),
  warn = TRUE,
  backend = c("auto", "dense", "banded"),
  cv.penalty.scale = c("training_mass", "none")
)

Arguments

Details

With training predictor locations normalized to ⁠[0,1]⁠, minimize

\frac{1}{2}\sum_i w_i(y_i-\eta_i)^2 +\lambda\|D^{(3)}f\|_1 +\frac{r_z}{2}\|\gamma\|^2 +\frac{r_f}{2\sum_i w_i}\sum_i w_i(f_i-\bar f_w)^2.

Here eta is the curve plus standardized covariates and offset. The intercept is unpenalized. Lambda is on the summed-loss scale; the default ridge penalties are zero. lambda=Inf restricts the curve to a global quadratic, not a constant. At lambda=0 with no covariates, offsets or ridge, predictions at distinct training locations are their weighted response means. Exact ties retain their individual loss contributions. Order is fixed at 2; no order-selection step is performed.

Gaussian refers to squared-error fitting. Fitting and CV do not estimate a noise variance, supply confidence intervals, or assume normality for the permutation assessment. Additional covariates and offsets are supported for fitting; assess.gaussian.trend.filtering supports only independent, unit-weight observations without adjustment.

Every CV split recomputes its training-only predictor range, knots and transformations. The default cv.penalty.scale="training_mass" scales lambda and ridge penalties by the training/full weight-sum ratio; "none" retains the earlier raw-penalty convention. This holds the mean-loss coefficients constant, not a common raw-coordinate operator across splits. Selection requires convergence on the full data and all folds. Equal CV scores select the first supplied lambda. Predictions outside the training range hold the curve at the nearest boundary; covariates and offsets continue to contribute. The banded backend scans all knots without thinning locations. Both solvers require finite arithmetic, stationarity, dual feasibility and L1 complementarity checks; failed fits remain explicit. Tolerances apply to mean weighted loss.

Value

A gaussian.trend.filtering.fit containing order=2L, the original data, penalty path, fits, convergence diagnostics, training design, and, when requested, cv.table with MSE and eligibility, fold-level loss sums and status, and selected.index. A fixed path has no selected index. cv.penalty.scale, cv.fold.lambda, cv.fold.penalty.factors, and cv.fold.ridge record the actual convention and training penalties. Standard prediction returns NA for unconverged candidates.

References

Wang, Smola and Tibshirani (2014). The Falling Factorial Basis and Its Statistical Applications. https://arxiv.org/abs/1405.0558

See Also

assess.gaussian.trend.filtering, fit.logistic.trend.filtering

Examples

x <- seq(0, 1, length.out = 40)
y <- sin(2*pi*x) + 0.1*cos(seq_along(x))
fit <- fit.gaussian.trend.filtering(x, y, lambda = c(Inf, 0.01),
    backend = "banded")
predict(fit, c(0.25, 0.75), lambda.index = 2)

Fit Graph Trend Filtering

Description

Fits graph trend filtering on a supplied undirected graph. Supported orders are 0L, 1L, and 2L:

\widehat\beta_\lambda = \arg\min_{\beta\in\mathbb R^n} \left\{ \frac{1}{2}\|y-\beta\|_2^2 + \lambda\|\Delta_w^{(k+1)}\beta\|_1 \right\}.

Usage

fit.graph.trend.filtering(
  adj.list,
  weight.list = NULL,
  y,
  order = 0L,
  lambda.grid = NULL,
  lambda.selection = c("cv", "fixed"),
  weight.rule = c("conductance", "sqrt.conductance", "unit"),
  operator.family = c("graph.laplacian.recursive", "path.divided.difference"),
  path.family = c("branch.continuation", "all.simple"),
  path.weighting = c("unit", "anchor.mean"),
  n.lambda = 40L,
  nfolds = 5L,
  foldid = NULL,
  maxsteps = 2000L,
  minlam = 0,
  approx = FALSE,
  rtol = 1e-07,
  btol = 1e-07,
  eps = 1e-04,
  verbose = FALSE,
  cv.penalty.scale = c("training_mass", "none")
)

Arguments

Details

Graph trend filtering is intended as an \ell_1-adaptive comparator to the \ell_2 graph smoothers now housed in geosmooth and the legacy rdgraph comparator. The low-pass smoother penalizes broad quadratic graph roughness through a term like

\eta f^\top L f,

which shrinks high-frequency variation everywhere. Graph trend filtering instead penalizes absolute graph differences. For order = 0L, this is

\lambda \sum_e \omega_e |\beta_j-\beta_i|,

so many edge differences can become exactly zero while selected edges carry jumps. In this sense order = 0L is locally adaptive: fitted values can be piecewise constant over graph regions. Higher orders use

\Delta_w^{(1)} = D_w,\qquad \Delta_w^{(2)} = L_w,\qquad \Delta_w^{(3)} = D_w L_w,

with L_w = D_w^\top D_w.

The supplied weight.list, when used, is interpreted as a graph trend-filtering edge weight or conductance for operator.family = "graph.laplacian.recursive". For operator.family = "path.divided.difference", weight.list is interpreted as metric edge length and is used to form local path coordinates. Use weight.rule = "sqrt.conductance" when the recursive graph square-root convention is desired, and weight.rule = "unit" for an unweighted recursive graph fused-lasso penalty.

The solver backend is currently genlasso. Cross-validation refits the generalized-lasso problem on vertex-held-out training sets using a selection matrix X; this is useful for small and moderate diagnostic problems but is not yet a scalable production graph-trend-filtering solver.

Value

A list of class "graph.trend.filtering.fit" containing fitted values, residuals, selected lambda, the operator, path metadata, and CV diagnostics when requested.

References

Wang, Y.-X., Sharpnack, J., Smola, A., and Tibshirani, R. J. (2016). Trend filtering on graphs. Journal of Machine Learning Research, 17(105), 1–41. https://www.jmlr.org/papers/v17/15-147.html

Examples

if (requireNamespace("genlasso", quietly = TRUE)) {
  adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
  fit.graph.trend.filtering(
    adj, y = c(0, 0.2, 1.8, 2), lambda.grid = 0.3,
    lambda.selection = "fixed", weight.rule = "unit"
  )
}

Fit a Local-Likelihood Density or Bernoulli Smoother

Description

Fits a local likelihood model at each evaluation point and reads off one raw fitted value at that evaluation anchor. Density and Bernoulli workflows should convert the returned fitted field with normalize.density when the target is a subject-occupation density.

Usage

fit.local.likelihood(
  X,
  y,
  X.eval = NULL,
  likelihood.family = c("density", "bernoulli"),
  support.size = min(15L, nrow(X)),
  degree = 1L,
  kernel = c("gaussian", "tricube", "epanechnikov", "triangular"),
  bandwidth.multiplier = 1,
  support.grid = NULL,
  degree.grid = NULL,
  kernel.grid = NULL,
  bandwidth.multiplier.grid = NULL,
  lambda.ridge.grid = NULL,
  foldid = NULL,
  cv.folds = 5L,
  cv.seed = 1L,
  coordinate.method = c("coordinates", "local.pca", "local.pls"),
  chart.dim = NULL,
  chart.dim.grid = NULL,
  selection.strategy = c("grid", "sparse_kd", "plateau_kd"),
  chart.dim.max = NULL,
  design.margin = 2L,
  auto.chart.support.metric = c("coordinates", "operator", "both"),
  auto.chart.selection.metric = c("coordinates", "operator"),
  quadrature.weights = NULL,
  lambda.ridge = 1e-08,
  min.local.mass = sqrt(.Machine$double.eps),
  min.nonzero.mass = 1L,
  fallback = c("degree0", "zero", "chart_kernel", "na"),
  optimizer = c("newton", "optim"),
  max.iter = 50L,
  tol = 1e-08,
  return.details = TRUE,
  chart.geometry = c("legacy", "euclidean", "linf.simplex"),
  unfold.faces = FALSE,
  intercept.mass = 0.5
)

Arguments

Details

The likelihood.family = "density" branch uses a local exponential tilt of the chart reference measure. The "bernoulli" branch uses a weighted local logistic likelihood and returns fitted probabilities at the evaluation anchors.

Value

A list with class "local_likelihood" containing fitted.values, selected controls, and local solver diagnostics.

Examples

X <- matrix(seq(0, 1, length.out = 21), ncol = 1)
y <- exp(-20 * (X[, 1] - 0.5)^2)
y <- y / sum(y)
fit.local.likelihood(X, y, likelihood.family = "density",
                     support.size = 7L, degree = 0L)

Fit Order-2 Logistic Trend Filtering

Description

Fits a piecewise-quadratic logit curve with an L1 penalty on its third divided differences. An R interface calls native C++ dense or banded solvers. No Python or external optimizer is required.

Usage

fit.logistic.trend.filtering(
  x,
  y,
  lambda,
  weights = NULL,
  covariates = NULL,
  offset = NULL,
  nuisance.ridge = 0,
  curve.ridge = 0,
  lambda.selection = c("fixed", "cv"),
  foldid = NULL,
  control = list(),
  warn = TRUE,
  backend = c("auto", "dense", "banded"),
  cv.penalty.scale = c("training_mass", "none")
)

Arguments

Details

Predictor locations are normalized to ⁠[0,1]⁠ using the training range. Lambda refers to third divided differences on that normalized scale, not on the original units of x. The minimized objective is

\sum_i w_i[\log(1+\exp(\eta_i))-y_i\eta_i] +\lambda\|D^{(3)}f\|_1 +\frac{r_z}{2}\|\gamma\|^2 +\frac{r_f}{2\sum_i w_i}\sum_i w_i(f_i-\bar f_w)^2.

Here eta is the curve plus standardized covariates and offset, gamma contains standardized covariate coefficients, and the curve includes an unpenalized intercept. Third differences annihilate quadratic polynomials. Both backends preserve the falling-factorial penalty exactly. The dense solver uses a scaled global basis and a logarithmic barrier. The banded solver uses locally supported linear combinations of the same functions, a fully corrective working set, banded curve Hessians, and a small Schur complement for covariates/centering. Every stopping check scans all knots; omitted working-set coefficients remain exactly zero and violating knots are admitted until the global KKT checks pass. It does not approximate the penalty or restrict the permitted knot locations. Exact ties share a curve value but retain their own likelihood contributions, covariates, and offsets. No jittering or approximate thinning is performed. Extra covariates can be confounded with the curve; without identifying constraints, their individual coefficients are not effects identifiable separately from the curve.

Convergence requires finite arithmetic, small stationarity residual, dual feasibility and small L1 complementarity. The reported dual vector lies in ⁠[-1,1]⁠; complementarity checks whether it is also an approximate subgradient at the returned coefficients. This is an approximate KKT check, not a separately computed primal-dual objective gap. Iteration limits and numerical failures remain explicit. Separation can prevent a finite minimizer, including in the unpenalized quadratic component. The classical Gaussian trend-filtering minimax result is not a guarantee for this model.

The dense backend is intended for moderate knot counts. The banded backend avoids observation-by-full-knot and full-knot-square arrays during fitting and prediction. Its storage is linear in observations and knots for a fixed number of covariates and penalties; computation also depends on working-set changes and Newton convergence. The conservative workspace estimate is checked before allocation. It is a solver workspace guard, not a process- wide memory limit. The quadratic initialization uses the dense solver only on the three polynomial columns and optional covariates.

CV recomputes training normalization, centering, and scaling in each fold. With the default cv.penalty.scale="training_mass", the penalty coefficients relative to weighted mean loss are held constant. The predictor range and standardized covariates still belong to each training split; scaling does not hold a common raw-coordinate operator fixed. A candidate with any failed training fit is ineligible for selection. All fits are refitted on the complete training data. Prediction outside its range uses the nearest boundary curve value; covariates and offsets continue to contribute normally. This boundary policy is not extrapolation of a quadratic polynomial.

Value

An object of class logistic.trend.filtering.fit, with fits, lambda, diagnostics, training transformations, optional cv.table and cv.fold.lambda, cv.fold.penalty.factors, cv.fold.ridge, the selected cv.penalty.scale, and selected.index. Failed fits are preserved for diagnosis but standard fitted/prediction methods return NA for them unless explicitly permitted.

See Also

assess.logistic.trend.filtering for a global constancy test and nested-CV predictive improvement for independent binary observations.

Examples

x <- seq(0, 1, length.out = 40)
y <- rep(c(0, 1, 0, 1, 1), 8)
model <- fit.logistic.trend.filtering(x, y, lambda = c(1, 0.1))
predict(model, seq(0, 1, length.out = 10), lambda.index = 2)

Fit A Local Polynomial Lifting Trend Filter

Description

Fits an LPL-TF estimator with a fixed operator:

\hat f = \arg\min_f \frac12\|y-f\|_2^2 + \lambda\|A_{\mathrm{LPL}}f\|_1.

Usage

fit.lpl.tf(
  X = NULL,
  y,
  operator = NULL,
  adj.list = NULL,
  weight.list = NULL,
  graph = NULL,
  degree = 2L,
  lambda.grid = NULL,
  lambda = NULL,
  lambda.selection = c("cv", "fixed"),
  operator.grid = NULL,
  foldid = NULL,
  cv.folds = 5L,
  cv.loss = c("mse", "rmse", "mae"),
  cv.seed = NULL,
  cv.repeats = 1L,
  solver = c("genlasso"),
  selection = c("min", "one.se"),
  n.lambda = 80L,
  maxsteps = 2000L,
  minlam = 0,
  approx = FALSE,
  rtol = 1e-07,
  btol = 1e-07,
  eps = 1e-04,
  verbose = FALSE,
  cv.penalty.scale = c("training_mass", "none"),
  ...
)

Arguments

Value

A list of class "lpl_tf".

Examples

if (requireNamespace("genlasso", quietly = TRUE)) {
  X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
  fit.lpl.tf(X, sin(2 * pi * X[, 1]), degree = 1L,
             support.type = "knn", support.size = 7L,
             lambda = 0.1, lambda.selection = "fixed")
}

Fit a Local Polynomial Smoother

Description

Fits a local polynomial smoother (LPS) and selects its support size, polynomial degree, and kernel by cross-validation. By default, the smoother works in the observed ambient coordinates: each prediction point uses its nearest training points in Euclidean distance, centers the support at the prediction point, fits a weighted local polynomial, and uses the fitted intercept as the prediction.

Usage

fit.lps(
  X,
  y,
  foldid = NULL,
  support.grid = c(10L, 15L, 20L),
  degree.grid = 0:2,
  kernel.grid = c("gaussian", "tricube"),
  cv.folds = 5L,
  cv.seed = 1L,
  X.eval = NULL,
  coordinate.method = c("coordinates", "local.pca"),
  chart.dim = NULL,
  chart.dim.grid = NULL,
  local.chart.method = c("pca", "second.order.svd"),
  auto.chart.support.metric = c("coordinates", "operator", "both"),
  auto.chart.selection.metric = c("coordinates", "operator"),
  backend = c("auto", "R", "cpp", "cpp.local.pca"),
  design.basis = c("orthogonal.polynomial.drop", "monomial", "weighted.qr",
    "weighted.qr.drop"),
  design.drop.tol = 1e-08,
  ridge.multiplier.grid = c(0, 1e-10, 1e-08),
  ridge.condition.max = 1e+12,
  unstable.action = c("na", "mean"),
  outcome.family = c("gaussian", "bernoulli", "binomial"),
  bandwidth.multiplier.grid = 1,
  keep.cv.predictions = FALSE,
  ridge.shrinkage.target = c("zero", "local.mean"),
  selection.strategy = c("grid", "sparse_kd", "plateau_kd"),
  chart.activation = c("none", "subject.od"),
  chart.activation.response = NULL,
  chart.activation.control = list(),
  chart.dim.max = NULL,
  design.margin = 2L
)

Arguments

Details

The optional coordinate.method = "local.pca" mode keeps the same support and kernel weighting rule, but builds the local polynomial in a local PCA chart centered at each prediction point. With chart.dim = "auto", the chart dimension is estimated as one global scalar from observed X only, using the same shared local-PCA dimension helper used by LPL-TF and S-LPL-TF. The experimental chart.dim = "local.auto" mode estimates an input-only local chart dimension separately for each prediction anchor.

Value

A list of class "lps" with response-scale fitted.values, unmodified local least-squares fitted.values.raw, selected parameters, a candidate CV table, the requested local.chart.method, and the effective chart method used for reporting. In outcome.family = "bernoulli" or "binomial" mode, fitted.values are response-scale probabilities in [0,1], fitted.values.raw are the un-clipped conditional-expectation estimates for "bernoulli" and the fitted probabilities for "binomial", cv.table$cv.brier.observed is the observed CV Brier score of the response-scale (clipped) predictions with Inf for candidates having any non-finite prediction, and probability.diagnostics records raw/clipped probability ranges, out-of-range fractions, and Brier/log-loss diagnostics. The Brier and log-loss diagnostics are defined only when the fitted predictions have the same length as the training response, which is the default X.eval = X path. In "binomial" mode, logistic.diagnostics records local logistic solve attempts, convergence statuses, fallback-path counts, event-rate fallback counts, and NA failure counts separately for CV and final fitting.

Choosing controls

First fit: supply X, y, foldid, and small support.grid, degree.grid, and kernel.grid candidates. See the installed function guide.

Selection: cv.folds, cv.seed, and selection.strategy control the search. Candidate combinations and folds multiply the fitting work; use lps.nested.cv() for an outer assessment of the tuning procedure.

Geometry: coordinate.method, chart.dim, and local chart controls describe neighborhoods and coordinate systems, not a different response model.

Numerics: backend, design.basis, rank thresholds and ridge controls govern stable local solves. Inspect lps.backend.diagnostics() before changing them. A numerical ridge is distinct from neighborhood selection.

Diagnostics: selected, cv.table, and diagnostics retain selection and local-fit information. A dense smoother.matrix needs one entry per evaluation/training pair; it is restricted to supported configurations.

References

Fan, J. and Gijbels, I. (1996). Local Polynomial Modelling and Its Applications. Chapman and Hall/CRC. ISBN 9780412983214.

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
y <- sin(2 * pi * X[, 1])
fit <- fit.lps(
  X, y, support.grid = 8L, degree.grid = 1L,
  kernel.grid = "tricube", cv.folds = 2L, backend = "R"
)
head(fit$fitted.values)

Fit Model-Averaged Local Polynomial Smoothing

Description

Fits the current model-averaged local polynomial smoother (MALPS) on a supplied coordinate matrix. MALPS fits local polynomial regressions around observed anchor points and averages all positive-weight local predictions at each training point. The current implementation supports supplied coordinates, observed or supplied coordinate anchors, fixed support parameters, deterministic cross-validation and exact dense GCV over support parameters, local PCA chart coordinates, ordinary and weighted new-response refits, coordinate-space new-point prediction, linear-smoother diagnostics, and fixed-selection bootstrap uncertainty summaries. Graph-geodesic supports are available for training and refit; graph-geodesic new-point prediction is deferred.

Usage

fit.malps(
  X,
  y,
  graph = NULL,
  adj.list = NULL,
  weight.list = NULL,
  graph.stage = "final",
  anchor.index = NULL,
  anchor.coordinates = NULL,
  degree = 2L,
  degree.grid = NULL,
  support.type = c("adaptive.radius", "knn", "fixed.radius"),
  support.size = NULL,
  support.grid = NULL,
  radius = NULL,
  radius.grid = NULL,
  min.support = NULL,
  min.support.grid = NULL,
  support.buffer = 3L,
  kernel = c("epanechnikov", "triangular", "gaussian", "tricube"),
  kernel.grid = NULL,
  model.weight.rule = c("none", "condition", "support", "boundary", "quality"),
  duplicate.action = c("keep", "error"),
  coordinate.method = c("coordinates", "local.pca"),
  chart.dim = NULL,
  support.metric = c("auto", "coordinates", "graph.geodesic"),
  auto.chart.support.metric = c("coordinates", "operator", "both"),
  auto.chart.selection.metric = c("coordinates", "operator"),
  support.selection = c("fixed", "cv", "gcv"),
  foldid = NULL,
  cv.folds = 5L,
  cv.loss = c("rmse", "mae", "mse"),
  cv.repeats = 1L,
  cv.seed = NULL,
  cv.one.se = FALSE,
  gcv.exact.max.n = 1000L,
  local.solver = c("auto", "normal.equations", "qr", "svd"),
  normal.equations.max.condition = 1e+08,
  robust.iterations = 0L,
  robust.tuning.constant = 6,
  verbose = FALSE,
  ...
)

Arguments

Details

Current development benchmarks suggest using the plain MALPS defaults for production-style fits unless a specific diagnostic motivates a non-default option: support.selection = "cv", robust.iterations = 0L, model.weight.rule = "none", and observed anchors (anchor.coordinates = NULL). Exact GCV selection (support.selection = "gcv") is a useful faster option for fixed-weight, non-robust MALPS fits, especially exploratory continuous smoothing runs, but it is not yet the universal default. Binary responses are treated as numeric conditional-expectation targets; fitted values are not clipped inside the smoother. Grid anchors, empirical-Bayes shrinkage, robust local residual reweighting, and conservative/smoothed GCV selection remain experimental controls rather than default recommendations.

Graph construction is supplied by dgraphs or by explicit adj.list/weight.list payloads. Coordinate MALPS paths and graph-geodesic payload validation are package-local; shortest-path distances are computed through dgraphs.

Value

A list of class "malps" with fitted values, local model coefficients, supports, averaging weights, diagnostics, and selection metadata.

Choosing controls

First fit: supply X, y, degree, and support.type with a suitable support.size or radius. Coordinate-distance fits support new predictions.

Selection: support and degree grids, foldid, and cv.folds govern held-out selection. Robust refits recompute residual weights; selection is not automatically repeated by refit.

Geometry: support.metric, graph inputs, anchors, and chart controls determine the local models. Graph-geodesic fits require aligned vertices and do not currently predict at new targets.

Numerics: local.solver, conditioning controls, and robust weighting have distinct meanings. Keep their defaults until diagnostics identify a need.

Diagnostics: inspect selected, model weights, and stored local fits. smoother.matrix allocates a dense matrix and applies only to fixed weighting; its max.n guard makes that cost explicit.

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
fit <- fit.malps(X, sin(2 * pi * X[, 1]), degree = 1L,
                 support.type = "knn", support.size = 8L)
head(fit$fitted.values)

Fit a Graph Field from Partially Observed Responses

Description

Estimate a conditional mean by weighted least squares with an unnormalized graph Laplacian penalty. Repeated observations at a vertex remain separate. Unobserved components return NA rather than an arbitrary constant.

Usage

fit.masked.graph.field(operator, vertex, y, lambda, weights = NULL)

Arguments

Value

A list with vertex fitted values, support flags, observation weights, objective, and relative linear-system residual.

Examples

op <- metric.graph.lowpass.operator(list(2L, c(1L, 3L), 2L),
                                   list(1, c(1, 1), 1))
fit <- fit.masked.graph.field(op, vertex = c(1, 3), y = c(2, 8), lambda = 0.5)
fit$fitted

Fit Metric-Conductance Graph Low-Pass Regression

Description

Fits graph-spectral low-pass regression using supplied conductances or transformed metric lengths. Missing-response heat fits use an observed-data penalized least-squares operator on a complete graph eigensystem.

Usage

fit.metric.graph.lowpass(
  adj.list,
  weight.list,
  y,
  conductance.rule = c("inverse.length.power", "exp.length", "exp.length.squared",
    "self.tuned.gaussian"),
  conductance.epsilon = 1e-08,
  conductance.alpha = 1,
  conductance.sigma = NULL,
  conductance.sigma.rule = c("edge.quantile", "median", "local.k"),
  conductance.sigma.quantile = 0.75,
  conductance.local.k = 5L,
  laplacian.type = c("unnormalized", "symmetric.normalized"),
  n.eigenpairs = 50L,
  filter.type = c("heat_kernel", "tikhonov", "cubic_spline", "gaussian", "exponential",
    "butterworth"),
  eta.grid = NULL,
  n.candidates = 40L,
  eigen.solver = c("auto", "sparse", "dense"),
  dense.eigen.threshold = 200L,
  dense.fallback.threshold = 5000L,
  dense.fallback = c("auto", "never", "always"),
  verbose = FALSE,
  eta.search = c("fixed", "guarded.gcv"),
  eta.expansion.factor = 3,
  eta.max.expansions = 3L,
  eta.identity.departure = 0.01,
  eta.truncation.tol = 1e-04,
  weight.type = c("length", "conductance"),
  spectrum = c("partial", "full"),
  numerical.policy = c("conservative", "verified.fallback")
)

Arguments

Details

Missing-response heat fits minimize observed squared error plus f^T (\exp(tL)-I) f. This reduces to ordinary heat filtering when all vertices are observed. GCV is q\,RSS/(q-\mathrm{tr}(S_O))^2, with q observed responses and observed-data operator S_O. Zero heat time is an identity diagnostic with undefined GCV and cannot be selected. A grid with no valid GCV candidate fails explicitly.

Value

A list of class "metric.graph.lowpass.fit". The GCV record includes per-candidate status, observed count and degrees of freedom. Missing-response residuals remain NA. The fitting.mode field distinguishes the observed heat-penalty fit from a fully observed spectral filter. For missing responses, the stored spectral filter weights describe the base heat kernel, not eigenvalues of the observed fitting matrix.

References

Gajer, P. and Ravel, J. (2025). Adaptive geometric regression for high-dimensional structured data. arXiv:2511.03817. doi:10.48550/arXiv.2511.03817

Examples

adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
lengths <- list(1, c(1, 1), c(1, 1), 1)
fit.metric.graph.lowpass(
  adj, lengths, y = c(0, 0.2, 1.8, 2), n.eigenpairs = 4L,
  eta.grid = c(0.1, 1), eigen.solver = "dense"
)

Nearest-Person and Nearest-Case Binary Probability Estimators

Description

Fit one of three estimators through a shared, single-threaded C++ backend. Ordinary nearest neighbors counts people; nearest-case estimators stop at the requested number of cases. Coordinates are processed one query at a time without constructing an all-pairs distance matrix.

Usage

fit.neighbor.probability(
  X = NULL,
  y,
  method = c("knn", "case_raw", "case_corrected"),
  k = NULL,
  k.grid = NULL,
  distances = NULL,
  metric = c("euclidean", "chebyshev", "linf.simplex"),
  normalization = c("none", "sum", "max"),
  id = NULL,
  group = NULL,
  foldid = NULL,
  selection = c("min", "one.se", "locality.one.se"),
  search = c("auto", "ann", "brute"),
  eps = 0
)

Arguments

Details

Rank answering training records by distance, then by the SHA-256 hash of paste0("nearest-cases-v1|", id). Hash collisions use ID order. Keep only the closest answering record per participant, before inspecting its outcome. Exclude the query's entire participant when its group is supplied.

If T_k is the rank including the kth case, the estimates are K^{-1}\sum_{j=1}^K y_{(j)}, k/T_k, and (k-1)/(T_k-1). No clipping, pseudocount or automatic baseline fallback is applied. Insufficient reachable people or cases return NA. The number of cases reachable can depend on the query when a participant's visits have discordant outcomes.

ANN expands candidate neighborhoods until the requested people or cases are available. Exact search completes the boundary before applying the ID tie rule. Approximate search sorts retrieved candidates by exact distance; for each retrieved participant it also checks all answering visits, keeping the closest. It can miss other participants and tied records. Approximate boundary diagnostics describe retrieved participants, not a complete shell. Available people/cases diagnostics still describe the full training set. Intrinsic surface search prepares sparse training profiles once per prediction call or CV split and streams a complete eligible distance row for each query. It uses dgraphs' shared intrinsic solver directly, without constructing a sample graph, limiting paths to endpoint faces, or truncating candidate neighbors. All three estimators share the same participant and tie rules. Per-pair worst-case work is quadratic in the number of features present in either profile. Distance workspace is linear in training records.

Trees are reused across queries and sizes, and across participants within each CV split. They are rebuilt for each prediction call and are not stored in the serializable fit. CV retains one tree or sparse intrinsic profile cache per split and repeat; a new prediction call rebuilds its cache. Native caches are released on return or error from prediction or tuning.

Tuning averages losses within participant, then across repeats, then across participants. Its standard error is the standard deviation of those participant losses divided by the square root of the number of people. Candidates with any unsupported validation prediction are unselectable. Numerical ties within 1e-14 favor larger sizes for minimum-error selection. The tuning score is not an independent performance estimate; use external participant-held outer folds for assessment.

This implementation matches the sPTB study's prediction arithmetic and minimum-error tuning when there is one record per participant, the same distances, IDs and folds. The repeated-record policy is an explicit extension. At a tied distance the prefix can split a tied shell, as in that study.

Value

A neighbor.probability.fit with selected k, method, backend = "C++", candidate scores, folds, and the metric, distance.backend, and training data needed for prediction. No training predictions are computed during fitting; predict(fit) computes leave-participant-out values.

Examples

X <- matrix(0:11, ncol = 1)
y <- rep(c(0, 1, 1), 4)
fit <- fit.neighbor.probability(X, y, method = "case_corrected", k = 2)
predict(fit, newdata = matrix(c(1.5, 5.5), ncol = 1))
tuned <- fit.neighbor.probability(X, y, method = "knn",
                                  foldid = rep(1:3, 4), k.grid = 1:4)
tuned$cv.table
# Intrinsic surface neighbors can differ from ambient Euclidean neighbors.
profiles <- rbind(c(0, 1, 0, 0), c(1, 1, 1, 1))
surface <- fit.neighbor.probability(profiles, c(1, 0), k = 1,
    metric = "linf.simplex", normalization = "max")
predict(surface, rbind(c(1, 0, 0, 0)))

Fit Prediction-Synchronized Local Polynomial Smoothing

Description

Experimental R reference implementation of prediction-synchronized local polynomial smoothing (PS-LPS). For fixed local polynomial parameters, PS-LPS fits one local polynomial chart per anchor and adds a quadratic penalty that synchronizes chart predictions on overlap points.

Usage

fit.ps.lps(
  X,
  y,
  foldid = NULL,
  support.size = NULL,
  degree = 2L,
  kernel = "gaussian",
  chart.dim = NULL,
  support.grid = NULL,
  degree.grid = NULL,
  kernel.grid = NULL,
  auto.chart.support.metric = c("coordinates", "operator", "both"),
  auto.chart.selection.metric = c("coordinates", "operator"),
  chart.dim.grid = NULL,
  selection.strategy = c("grid", "sparse_kd", "plateau_kd"),
  chart.dim.max = NULL,
  design.margin = 2L,
  lambda.sync.grid = c(0, 0.001, 0.01, 0.1, 1, 10),
  lambda.sync.search = c("grid", "guarded"),
  lambda.sync.selection = c("cv", "fixed"),
  local.candidate.search = c("screened", "full", "subgrid"),
  local.candidate.search.control = list(),
  lambda.sync.search.control = list(),
  lambda.diagnostics = c("all", "selected"),
  lambda.ridge = 1e-08,
  design.basis = c("monomial", "weighted.qr", "weighted.qr.drop",
    "orthogonal.polynomial.drop"),
  design.drop.tol = sqrt(.Machine$double.eps),
  ridge.multiplier.grid = NULL,
  ridge.condition.max = Inf,
  sync.neighbor.size = NULL,
  overlap.weight = c("normalized.product", "product"),
  chart.activation = c("none", "subject.od"),
  chart.activation.response = NULL,
  chart.activation.control = list(),
  ps.lps.geometry.cache = NULL,
  ps.lps.local.pca.supports = NULL,
  cv.folds = 5L,
  cv.seed = 1L
)

Arguments

Value

A list with fitted values, selected lambda, CV table, diagnostics, and fitted chart coefficients.

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
fit.ps.lps(
  X, sin(2 * pi * X[, 1]), support.size = 7L, degree = 1L,
  chart.dim = 1L, lambda.sync.grid = 0,
  lambda.sync.selection = "fixed", cv.folds = 2L
)

Fit A Parallel-Transport Trend Filter From A PTTF Operator

Description

Fits experimental PTTF regression models from a transported difference operator produced by pttf.operator. The implementation focuses on one response vector and explicit operator-row policies, especially the 1-D boundary policy used to separate interior transported rows from endpoint rows.

Usage

fit.pttf.trend.filtering(
  geometry = NULL,
  operator = NULL,
  X = NULL,
  y,
  derivative.order = 3L,
  penalty = c("l1", "l2"),
  lambda.grid = NULL,
  lambda.selection = c("cv", "fixed"),
  n.lambda = 80L,
  nfolds = 5L,
  foldid = NULL,
  cv.loss = c("mse", "mae"),
  selection = c("min", "one.se"),
  lambda.extension = c("none", "auto"),
  solver = c("genlasso", "admm", "auto"),
  operator.row.policy = c("all", "drop.line.boundary", "diagnostic.only"),
  line.order = NULL,
  boundary.trim = NULL,
  row.mass.rule = c("none", "node.mass"),
  row.normalize = c("none", "l2"),
  weights = NULL,
  maxsteps = 2000L,
  minlam = 0,
  approx = FALSE,
  rtol = 1e-07,
  btol = 1e-07,
  eps = 1e-04,
  admm.rho = NULL,
  admm.maxiter = 2000L,
  admm.abstol = 1e-04,
  admm.reltol = 0.001,
  verbose = FALSE,
  admm.adaptive.rho = TRUE,
  ...
)

Arguments

Value

A list of class "pttf.trend.filtering.fit".

Examples

n <- 10L
X <- matrix(seq(0, 1, length.out = n), ncol = 1)
adj <- lapply(seq_len(n), function(i) intersect(c(i - 1L, i + 1L), 1:n))
lengths <- Map(function(i, j) abs(X[j, 1] - X[i, 1]), seq_len(n), adj)
geometry <- pttf.geometry(
  X, adj, lengths, graph = "supplied", tangent.dim = 1L
)
operator <- pttf.operator(geometry, derivative.order = 2L)
fit.pttf.trend.filtering(
  operator = operator, y = sin(2 * pi * X[, 1]), penalty = "l2",
  lambda.grid = 0.1, lambda.selection = "fixed"
)

Fit A Fixed-Operator Synchronized LPL-TF Model

Description

Fits the fixed-operator S-LPL-TF estimator

\hat f = \arg\min_f \frac12\|y-f\|_2^2 + \lambda_1\|A_{\mathrm{LPL}}f\|_1 + \frac{\lambda_2}{2}\|C_{\mathrm{sync}}f\|_2^2.

Usage

fit.slpl.tf(
  X = NULL,
  y,
  operator = NULL,
  lambda1 = NULL,
  lambda2 = 0,
  lambda1.grid = NULL,
  lambda2.grid = NULL,
  lambda.selection = c("fixed", "cv"),
  operator.grid = NULL,
  foldid = NULL,
  cv.folds = 5L,
  cv.loss = c("mse", "rmse", "mae"),
  cv.seed = NULL,
  cv.repeats = 1L,
  solver = c("genlasso"),
  selection = c("min", "one.se"),
  maxsteps = 2000L,
  minlam = 0,
  approx = FALSE,
  rtol = 1e-07,
  btol = 1e-07,
  eps = 1e-04,
  verbose = FALSE,
  ...
)

Arguments

Value

A list of class "slpl_tf".

Examples

if (requireNamespace("genlasso", quietly = TRUE)) {
  X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
  fit.slpl.tf(X, sin(2 * pi * X[, 1]), degree = 1L,
              support.type = "knn", support.size = 7L,
              lambda1 = 0.1, lambda2 = 0, lambda.selection = "fixed")
}

Fit SSRHE-Style Hessian L1 Regression

Description

Fits an \ell_1-adaptive SSRHE Hessian comparator using the row operator A from ssrhe.hessian.operator:

\widehat f_\lambda = \arg\min_f \left\{ \frac{1}{2}\sum_i w_i(y_i-f_i)^2 + \lambda \|Af\|_1 \right\}.

Usage

fit.ssrhe.hessian.l1.regression(
  X,
  y,
  k = NULL,
  tangent.dim,
  lambda.grid = NULL,
  lambda.selection = c("cv", "fixed"),
  weights = NULL,
  n.lambda = 40L,
  nfolds = 5L,
  fold.id = NULL,
  loss = c("mse", "mae"),
  selection = c("min", "one.se"),
  nn.index = NULL,
  neighborhood.type = c("knn", "adaptive.radius", "supplied"),
  support.index = NULL,
  adaptive.k.scale = NULL,
  radius.rule = c("geomean", "max", "min"),
  radius.factor = 1.25,
  min.support = NULL,
  max.support = NULL,
  support.buffer = 2L,
  support.topup = c("nearest", "none"),
  tangent.dim.rule = c("fixed", "eigen.cumulative"),
  eigen.tolerance = 0.95,
  derivative.order = 2L,
  pinv.tol = sqrt(.Machine$double.eps),
  local.solver = c("auto", "normal.equations", "svd", "qr"),
  normal.equations.max.condition = 10000,
  solver = c("genlasso", "admm", "auto"),
  row.scaling = c("none", "l2"),
  admm.rho = NULL,
  admm.maxiter = 2000L,
  admm.abstol = 1e-06,
  admm.reltol = 1e-05,
  maxsteps = 2000L,
  minlam = 0,
  approx = FALSE,
  rtol = 1e-07,
  btol = 1e-07,
  eps = 1e-04,
  support.selection = c("rule", "cv"),
  support.grid = NULL,
  support.cv.max.candidates = 8L,
  return.local.diagnostics = FALSE,
  return.timing = FALSE,
  verbose = FALSE,
  admm.adaptive.rho = TRUE
)

Arguments

Details

This function exposes the SSRHE local polynomial derivative estimator as an \ell_1 penalty, rather than as the original SSRHE \ell_2 Hessian-energy penalty used by fit.ssrhe.hessian.regression. With derivative.order = 2L, A contains local quadratic/Hessian rows. With derivative.order = 3L, A contains local cubic third-derivative tensor rows with the same symmetric component scaling used by ssrhe.hessian.operator. The existing \ell_2 fit solves a linear system involving B=A^\top A. This function instead keeps the rows of A and penalizes their absolute values. It is therefore closer in spirit to graph trend filtering, while retaining SSRHE's local-PCA derivative construction.

The default solver backend is genlasso. For larger or numerically fragile third-order operators, solver = "admm" fits the requested lambda values directly without computing a generalized-lasso path, and all backends use ADMM for the returned coefficients. solver = "auto" can additionally recover from path construction errors with an explicit grid. Cross-validation removes held-out labels from the data-fit term by setting their weights to zero, then scores predictions on those held-out labels.

ADMM uses proximal numerical stabilization, which does not add a ridge penalty to the requested objective. Convergence requires primal feasibility, the dual stopping residual, and stationarity of the original weighted L1 objective to meet their reported tolerances. A finite fixed-lambda iterate that reaches its cap is returned with a warning and solver$converged = FALSE; it is not a completed fit. Inspect solver$admm before using it, and consider a larger iteration limit or an adaptive penalty. Non-finite selected fits are errors.

CV candidates are eligible only when every fold has an acceptable solve. ADMM iterates that do not converge are not scored. The CV result retains fold.status, fold.converged, fold.diagnostics, n.successful.folds, and eligible. If no candidate qualifies, an error of class ssrhe_l1_cv_error retains fold diagnostics for inspection with tryCatch. A selected full-data ADMM fit that fails to converge is also an error in CV mode.

A generalized-lasso path alone can be inaccurate for nearly rank-deficient operators. It is retained in path for diagnostics, but its coefficients are not the returned fit. solver$backend = "genlasso_admm_refined" identifies this route; solver$path.objective records the unrefined objectives at the grid values. Backend warnings, including its numerical ridge for underdetermined designs, are retained in solver$warnings. CV uses the original objective without that additional ridge. The default ADMM tolerances are absolute 1e-6 and relative 1e-5; convergence describes numerical stopping at those tolerances, not predictive accuracy. With row.scaling = "l2", the reported objective uses the scaled penalty; energies$penalty.l1 is that penalty norm, while energies$hessian.l1 retains the unscaled operator norm.

With support.selection = "cv", each adaptive-radius support candidate builds a fresh SSRHE operator and runs the usual lambda CV with shared fold assignments. The selected fit is the candidate with the smallest selected CV error. This is useful when the default support-size rule is too rigid, but it can be much slower than tuning lambda for one fixed operator.

Value

A list of class "ssrhe.hessian.l1.fit" containing fitted values, residuals, selected lambda, the SSRHE operator, generalized-lasso path metadata, lambda-grid fitted values, and CV diagnostics when requested.

Choosing controls

First fit: specify X, y, k, tangent.dim, and lambda.grid. For a single supplied penalty set lambda.selection = "fixed" and solver = "admm" to avoid constructing a path.

Selection: fold.id, nfolds, loss, and selection control label CV. A grid requires one solve per penalty per fold; support CV also rebuilds geometry. Incomplete candidates are excluded, not silently scored.

Geometry: neighborhood, tangent and derivative controls define A. row.scaling changes the penalty itself and thus the statistical objective.

Numerics: ⁠admm.*⁠ controls stopping and algorithmic penalties for all backends. Path controls (maxsteps, rtol, btol, eps) affect the optional proposal path. They do not certify the final original-objective solve.

Diagnostics: solver$admm, candidate/fold status and objective components describe numerical acceptance. path and beta.grid can use substantial memory; use a compact explicit grid with ADMM when no full path is needed.

References

Kim, K. I., Steinke, F., and Hein, M. (2009). Semi-supervised regression using Hessian energy with an application to semi-supervised dimensionality reduction. Advances in Neural Information Processing Systems 22. https://papers.nips.cc/paper_files/paper/2009/hash/f4552671f8909587cf485ea990207f3b-Abstract.html

Examples

X <- matrix(seq(0, 1, length.out = 12), ncol = 1)
fit.ssrhe.hessian.l1.regression(
  X, sin(2 * pi * X[, 1]), k = 6L, tangent.dim = 1L,
  lambda.grid = 0.05, lambda.selection = "fixed", solver = "admm"
)

Fit SSRHE-Style Hessian-Energy Regression

Description

Fits the \ell_2 Hessian-energy regularized estimator associated with ssrhe.hessian.operator. This is a direct SSRHE-style comparator for Hessian smoothing on point clouds; it is not a replacement for graph trend-filtering regression and is distinct from \ell_1-adaptive graph trend filtering.

Usage

fit.ssrhe.hessian.regression(
  X,
  y,
  k = NULL,
  tangent.dim,
  lambda1,
  lambda2 = 0,
  weights = NULL,
  nn.index = NULL,
  neighborhood.type = c("knn", "adaptive.radius", "supplied"),
  support.index = NULL,
  adaptive.k.scale = NULL,
  radius.rule = c("geomean", "max", "min"),
  radius.factor = 1.25,
  min.support = NULL,
  max.support = NULL,
  support.buffer = 2L,
  support.topup = c("nearest", "none"),
  tangent.dim.rule = c("fixed", "eigen.cumulative"),
  eigen.tolerance = 0.95,
  derivative.order = 2L,
  stabilizer = lambda2 > 0,
  pinv.tol = sqrt(.Machine$double.eps),
  local.solver = c("auto", "normal.equations", "svd", "qr"),
  normal.equations.max.condition = 10000,
  ridge = 0,
  return.A = TRUE,
  return.local.diagnostics = FALSE,
  return.timing = FALSE,
  verbose = FALSE,
  lambda.selection = c("fixed", "cv", "gcv"),
  lambda1.grid = NULL,
  lambda2.grid = NULL,
  cv.control = list(),
  gcv.control = list()
)

Arguments

Details

For each response column, this function solves

(W + \lambda_1 B + \lambda_2 B_S + \epsilon I)\hat f = Wy,

where W is the diagonal matrix of observation weights and \epsilon is ridge. Fully observed data with lambda1 = lambda2 = ridge = 0 therefore reproduce the observed response exactly. Missing responses are excluded from the data-fit term by setting their weights to zero. The Matlab SSRHE semi-supervised convention is represented by a 0/1 labeled-indicator weights vector, or equivalently by setting unlabeled responses to NA: labeled vertices contribute unit diagonal data-fit terms and unlabeled vertices contribute no data-fit term.

The lower-level operator is matched to the Kim–Steinke–Hein SSRHE Matlab construction: self-including kNN neighborhoods, local PCA charts, a fixed-intercept local quadratic fit, doubled diagonal Hessian components, and optional supplemental stabilizer. The package exposes both A and B=A^\top A; the fitted \ell_2 estimator uses B.

Value

A list of class "ssrhe.hessian.fit" containing fitted values, residuals, input response, weights, lambda parameters, objective/energy diagnostics, the reused operator, solver metadata, the call, and lambda.selection. CV adds class "ssrhe.hessian.cv.fit", cv.table, fold.id, and selection; GCV adds class "ssrhe.hessian.gcv.fit", gcv.table, gcv.trace, and selection. Support searches retain selected.support and the corresponding support.cv.table or support.gcv.table.

Penalty selection

Fixed fitting supports response vectors or matrices. CV and GCV currently select penalties for one response vector at a time. Label CV removes each validation fold from the data-fit term and scores its held-out observed positive-weight labels. It supports missing responses. GCV requires a fully observed response and strictly positive observation weights, and minimizes

\frac{n^{-1}\sum_i w_i(y_i-\hat f_i)^2} {(1-\mathrm{tr}(S_\lambda)/n)^2},

where S_\lambda=(W+\lambda_1 B+\lambda_2 B_S+\epsilon I)^{-1}W. The exact trace is deterministic. The Hutchinson option uses Rademacher probes, reports trace standard errors, and is reproducible with gcv.control = list(trace.method = "hutchinson", trace.seed = ...).

Both selection modes can also search adaptive-radius support profiles: set support.selection to the same mode in its control list. Each row of support.grid contains adaptive.k.scale, min.support, and optionally max.support. A NULL grid is constructed by ssrhe.support.grid. Each profile requires a new operator and a full penalty search, so this can substantially increase runtime.

Migration

The former fit.ssrhe.hessian.regression.cv() and fit.ssrhe.hessian.regression.gcv() entry points have been retired. Use this function with lambda.selection = "cv" or "gcv". Move fold.id to cv.control$foldid, nfolds to cv.control$cv.folds, and loss/selection to cv.control. Move gcv.trace.method, gcv.trace.n.probes, and gcv.trace.seed to trace.method, trace.n.probes, and trace.seed in gcv.control. Support search settings also belong in the active control list; both former support candidate limits are now named support.max.candidates.

Choosing controls

First fit: supply X, y, neighborhood size k, tangent.dim, and lambda1; use lambda.selection = "fixed" for a specified penalty.

Selection: lambda.selection, lambda1.grid, and lambda2.grid select penalties. Fold settings belong in cv.control; trace settings belong in gcv.control. Label CV and GCV are different criteria. Support selection rebuilds the operator for each candidate and multiplies this work.

Geometry: neighborhood, adaptive-radius, tangent and derivative controls determine the operator. Reuse it when fitting multiple responses.

Numerics: ridge and local conditioning controls stabilize solves. Exact GCV traces require more solves than stochastic trace estimates; the latter report sampling uncertainty and need a reproducible probe seed.

Diagnostics: inspect the selection table, energies, observed-label mask, and solver metadata. fitted() and residuals() preserve response matrices.

References

Kim, K. I., Steinke, F., and Hein, M. (2009). Semi-supervised regression using Hessian energy with an application to semi-supervised dimensionality reduction. Advances in Neural Information Processing Systems 22. https://papers.nips.cc/paper_files/paper/2009/hash/f4552671f8909587cf485ea990207f3b-Abstract.html

Examples

X <- matrix(seq(0, 1, length.out = 12), ncol = 1)
fit.ssrhe.hessian.regression(
  X, sin(2 * pi * X[, 1]), k = 6L, tangent.dim = 1L, lambda1 = 0.1
)
fit.ssrhe.hessian.regression(
  X, sin(2 * pi * X[, 1]), k = 6L, tangent.dim = 1L,
  lambda.selection = "cv", lambda1.grid = c(0.01, 0.1),
  cv.control = list(cv.folds = 2L)
)
fit.ssrhe.hessian.regression(
  X, sin(2 * pi * X[, 1]), k = 6L, tangent.dim = 1L,
  lambda.selection = "gcv", lambda1.grid = c(0.01, 0.1)
)

Fit Subject-Occupation Density

Description

Convenience wrapper that converts subject visit indices into a density input. Density-native methods dispatch to fit.density. LPS/PS-LPS methods fit ordinary smoother objects first and then call normalize.density; they are subject-occupation workflows, not standalone density-native methods.

Usage

fit.subject.od(
  X,
  subject.index,
  method = c("empirical", "graph_random_walk", "lps_count", "ps_lps_count",
    "lps_logistic_binary", "chart_kernel", "local_likelihood_density",
    "local_likelihood_bernoulli"),
  graph = NULL,
  graph.control = list(),
  od.control = list(),
  return.details = TRUE,
  od.cv = c("none", "visit"),
  visit.foldid = NULL,
  visit.cv.folds = 5L,
  visit.cv.seed = 1L,
  visit.cv.epsilon = 1e-15,
  ...
)

Arguments

Details

The "lps_logistic_binary" method name is historical. In the current OD workflow it calls fit.lps(..., outcome.family = "bernoulli"), which fits a clipped probability field with the LPS least-squares core and then converts that field to a density with normalize.density. It is not the outcome.family = "binomial" local-logistic IRLS path.

When od.cv = "visit", graph random-walk candidates may pass OD-level grids inside graph.control: walk.step.grid, affinity.method.grid, affinity.scale.grid, affinity.epsilon.grid, and normalize.grid. Missing affinity.scale values mean the usual data-derived affinity scale.

When od.cv = "visit", chart-based and LPS-family methods may pass chart.dim.grid through .... For "chart_kernel", "local_likelihood_density", "local_likelihood_bernoulli", "lps_count", "lps_logistic_binary", and "ps_lps_count", this grid compares fixed integer chart dimensions with optional "auto" and "local.auto" chart-dimension policies under the same held-out-visit negative-log-occupation score. For LPS-family OD visit CV, each outer candidate is passed to the source smoother as a scalar local model configuration; this avoids nested row-level multi-candidate selection inside each held-out-visit fold.

Value

A list of class "density_fit". Its rho component is the estimated probability mass over the rows of X, and empirical.rho is the normalized visit-count mass. The subject component summarizes the number of visits, the number of distinct visited support points, the largest visit multiplicity, and the fraction of repeatedly visited support points. Other components describe the fitted method, normalization accounting, smoothing diagnostics, and warnings as documented in fit.density. With od.cv = "visit" and return.details = TRUE, the result also contains the visit-level cross-validation table, fold assignments, and held-out predicted masses.

Examples

X <- matrix(seq(0, 1, length.out = 6), ncol = 1)
fit.subject.od(X, subject.index = c(1L, 3L, 3L, 6L),
               method = "empirical")

Get Boundary Vertices of a Region in a Graph

Description

Identifies the boundary vertices of a specified region in a graph, defined as vertices in the region that have at least one neighbor outside the region.

Usage

get.region.boundary(adj.list, region)

Arguments

Details

This function determines which vertices in a specified region are positioned at the boundary - meaning they have at least one neighbor that is not part of the region. These boundary vertices are often treated differently in graph-based algorithms, such as harmonic smoothing, where their values are typically fixed as constraints.

The boundary definition used matches the one implemented in the C++ harmonic smoothing functions, focusing on vertices that have connections to the outside of the region.

Value

An integer vector containing the indices of the boundary vertices, which is a subset of the input region vector. Returns an empty integer vector if no boundary vertices exist.

See Also

perform.harmonic.smoothing, harmonic.smoother

Examples

# Create a simple grid graph adjacency list (5x5 grid)
create_grid_adj_list <- function(n_rows, n_cols) {
  n <- n_rows * n_cols
  adj_list <- vector("list", n)
  for (i in 1:n_rows) {
    for (j in 1:n_cols) {
      v <- (i-1) * n_cols + j
      neighbors <- numeric(0)

      # Add neighbors (up, down, left, right)
      if (i > 1) neighbors <- c(neighbors, (i-2) * n_cols + j)  # up
      if (i < n_rows) neighbors <- c(neighbors, i * n_cols + j)  # down
      if (j > 1) neighbors <- c(neighbors, (i-1) * n_cols + (j-1))  # left
      if (j < n_cols) neighbors <- c(neighbors, (i-1) * n_cols + (j+1))  # right

      adj_list[[v]] <- neighbors
    }
  }
  return(adj_list)
}

# Create a 5x5 grid graph
grid_adj_list <- create_grid_adj_list(5, 5)

# Define a region (central 3x3 subgrid)
central_region <- c(7:9, 12:14, 17:19)

# Find boundary vertices of the central region
boundary <- get.region.boundary(grid_adj_list, central_region)
print(boundary)
# Expected output: c(7, 8, 9, 12, 14, 17, 18, 19)
# (all except the center vertex 13)

Construct a Graph Trend-Filtering Operator

Description

Builds the weighted graph difference operators used by graph trend filtering. The default operator.family = "graph.laplacian.recursive" uses recursive graph operators

k=0:\quad \Delta_w^{(1)} = D_w,

k=1:\quad \Delta_w^{(2)} = L_w = D_w^\top D_w,

and

k=2:\quad \Delta_w^{(3)} = D_w L_w.

The experimental operator.family = "path.divided.difference" builds canonical path rows and computes position-aware finite differences along graph paths using weight.list as metric edge lengths.

Usage

graph.trend.filtering.operator(
  adj.list,
  weight.list = NULL,
  order = 0L,
  weight.rule = c("conductance", "sqrt.conductance", "unit"),
  operator.family = c("graph.laplacian.recursive", "path.divided.difference"),
  path.family = c("branch.continuation", "all.simple"),
  path.weighting = c("unit", "anchor.mean"),
  return.sparse = TRUE
)

Arguments

Value

A list of class "graph.trend.filtering.operator" containing the canonical graph, edge table, trend-filtering weights, incidence, weighted Laplacian, selected penalty operator, and nullity estimate.

Examples

adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
graph.trend.filtering.operator(adj, order = 1L, weight.rule = "unit")

Perform Harmonic Smoothing with Topology Tracking

Description

Relaxes graph values while recording changes in local extrema. Boundary vertices have neighbors outside the region and retain their original values.

Usage

harmonic.smoother(
  adj.list,
  weight.list,
  values,
  region.vertices,
  max.iterations = 100,
  tolerance = 1e-06,
  record.frequency = 1,
  stability.window = 3,
  stability.threshold = 0.05
)

Arguments

Details

Uses the same inverse-length neighbor averages and boundary definition as perform.harmonic.smoothing. Unlike the basic function, it also relaxes regions with no boundary. Such a problem is not uniquely anchored; synchronous updates can oscillate on bipartite graphs, including paths.

The recorded diagnostic compares sets of (vertex, minimum/maximum) pairs. Its value is one minus twice the number of shared pairs divided by the sum of the two set sizes (zero for two empty sets). A detection means that the last stability.window recorded differences are at most stability.threshold. It does not establish numerical convergence, optimal smoothing, noise removal, or stability of full attraction basins. In particular, recording every two iterations can conceal an alternating oscillation. Always inspect converged and max_residual separately. Stability detection does not stop relaxation or select the returned fit. The initial and final states are always recorded, including a final partial recording interval. Disconnected components without fixed boundary values have the same lack of anchoring as a region without a boundary.

Value

A list of class "harmonic_smoother" containing:

See Also

perform.harmonic.smoothing for basic smoothing, plot.harmonic_smoother for visualization methods, summary.harmonic_smoother for summary statistics

Examples

adj.list <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
weight.list <- list(1, c(1, 1), c(1, 1), 1)
values <- c(0, 2, -1, 1)
result <- harmonic.smoother(
  adj.list, weight.list, values, region.vertices = seq_along(values),
  max.iterations = 20, record.frequency = 2
)
smoothed.values <- result$harmonic_predictions
smoothed.values

Construct a Local Polynomial Lifting Trend-Filtering Operator

Description

Builds the Local Polynomial Lifting Trend Filtering (LPL-TF) analysis operator from local polynomial prediction residuals. For each observed target point i, the operator predicts f_i from nearby values f_j using a local polynomial model and stores the residual row

r_i(f) = f_i - \sum_{j \in S_i} h_{ij} f_j.

Usage

lpl.tf.operator(
  X,
  adj.list = NULL,
  weight.list = NULL,
  graph = NULL,
  graph.stage = "final",
  anchor.index = NULL,
  anchor.coordinates = NULL,
  degree = 2L,
  support.type = c("adaptive.radius", "knn", "fixed.radius"),
  support.size = NULL,
  radius = NULL,
  min.support = NULL,
  support.buffer = 3L,
  kernel = c("epanechnikov", "triangular", "gaussian", "tricube"),
  coordinate.method = c("coordinates", "local.pca"),
  chart.dim = NULL,
  support.metric = c("auto", "coordinates", "graph.geodesic"),
  auto.chart.support.metric = c("coordinates", "operator", "both"),
  auto.chart.selection.metric = c("coordinates", "operator"),
  exclude.self = TRUE,
  row.normalize = c("l2", "none", "l1"),
  local.solver = c("auto", "normal.equations", "qr", "svd"),
  normal.equations.max.condition = 1e+08,
  duplicate.action = c("keep", "error"),
  drop.rank.deficient = TRUE,
  verbose = FALSE,
  ...
)

Arguments

Details

The operator uses observed anchors, excludes the target point from every predictor support, and drops requested-degree rank-deficient rows with explicit metadata. There is no silent degree downgrade. Graph support distances affect support selection and kernel bandwidths only; local polynomial coordinates are controlled separately by coordinate.method. Graph construction is supplied by dgraphs or by explicit adj.list/weight.list payloads. Coordinate support paths and graph-geodesic payload validation are package-local; shortest-path distances are computed through dgraphs.

If a_i is the complete assembled residual row, then row.normalize = "l2" uses a_i / \|a_i\|_2, and row.normalize = "l1" uses a_i / \|a_i\|_1. Normalization is applied to the whole row, not only to predictor coefficients.

Value

A list of class "lpl_tf_operator" containing the sparse operator matrix A, row metadata, supports, local design summaries, diagnostics, settings, and the matched call.

Examples

X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
lpl.tf.operator(X, degree = 1L, support.type = "knn",
                support.size = 7L, kernel = "gaussian")

Report LPS Backend and Chart-Dimension Diagnostics

Description

Builds a compact one-row diagnostic table for a fitted local polynomial smoother. The table records the requested backend, the backend actually used, the requested chart-dimension rule, the resolved chart dimension, selected tuning parameters, and whether the fit follows the current deployable local-PCA auto-dimension contract.

Usage

lps.backend.diagnostics(object)

Arguments

Details

The current backend policy is conservative: backend = "auto" uses the C++ backend for ambient-coordinate LPS, but uses the R reference backend for coordinate.method = "local.pca". The native local-PCA backend "cpp.local.pca" remains an explicit opt-in backend. This helper makes that policy visible in reports and downstream experiment manifests without changing the default.

For real-data local-PCA runs, the deployable chart-dimension contract is the ordinary local.chart.method = "pca" path with chart.dim = "auto" or "local.auto" and observed-covariate auto-dimension diagnostics. In P7-style experiments this is paired with auto.chart.support.metric = "both" and auto.chart.selection.metric = "operator". The experimental "second.order.svd" chart path is reported explicitly but is not certified by this deployable local-PCA contract. For LPS itself, which uses coordinate supports, the operator-support diagnostic is currently equivalent to the coordinate-support diagnostic; the fields are still recorded so the same manifest schema can be shared with LPL-TF and S-LPL-TF experiments.

Value

A one-row data.frame with backend, chart-dimension, selection, candidate-count, and policy fields.

Examples

X <- matrix(seq(0, 1, length.out = 16), ncol = 1)
fit <- fit.lps(X, X[, 1]^2, support.grid = 6L, degree.grid = 1L,
               kernel.grid = "gaussian", cv.folds = 2L, backend = "R")
lps.backend.diagnostics(fit)

Build a Cluster-Respecting Fold Assignment

Description

Assigns whole clusters to cross-validation folds so that no cluster is split across folds: every observation of a cluster receives the same fold id, and therefore no training set ever shares a cluster with the fold held out against it. Folds are built by a deterministic size-balanced greedy rule: clusters are taken largest first (ties by first appearance in cluster.id) and each is placed into the currently smallest fold (ties by lowest fold index). With v equal to the number of clusters this reduces to leave-cluster-out.

Usage

lps.grouped.foldid(cluster.id, v = 5L, shuffle.seed = NULL)

Arguments

Details

The assignment is fully deterministic when shuffle.seed is NULL. Supplying shuffle.seed applies one seeded permutation to the cluster order before the greedy pass (the seed is consumed immediately via set.seed), giving a reproducible randomized variant; callers following the LPS evidence conventions should record the seed they pass.

Value

An integer fold-id vector of length length(cluster.id) with values in 1:v; every fold is nonempty and every cluster maps to exactly one fold.

See Also

lps.nested.cv() for nested cross-validation that can consume grouped folds at both the outer and inner level.

Examples

cluster <- rep(letters[1:6], each = 2)
lps.grouped.foldid(cluster, v = 3L)

Nested Cross-Validation for LPS with Explicit, Recorded Folds

Description

Runs outer-fold nested cross-validation around fit.lps(): for each outer fold, candidate selection is one ordinary fit.lps call on the inner-training rows only, with an explicit inner fold id and X.eval set to the held-out outer-test rows, so the held-out fold never participates in inner selection. The pooled outer-test error is the nested generalization estimate. The same call also computes the selected-min arm — an ordinary fit.lps on all rows using the same outer.foldid — so optimism comparisons are paired on one fold assignment by construction.

Usage

lps.nested.cv(
  X,
  y,
  outer.foldid,
  fit.args = list(),
  inner.folds = 5L,
  cluster.id = NULL,
  inner.foldid.method = c("round.robin", "grouped"),
  inner.shuffle.seed = NULL
)

Arguments

Details

The existing fit.lps behavior is untouched: this utility consumes the public API with explicit foldid only. It currently supports outcome.family = "gaussian" (the default) and scores with RMSE.

Inner folds are built per outer fold and fully recorded in the return value: "round.robin" assigns rep_len(1:inner.folds, n) over the inner-training rows in row order (deterministic; if inner.shuffle.seed is supplied, fold k — by position — uses set.seed(inner.shuffle.seed + k) and permutes the balanced assignment, giving a reproducible randomized variant); "grouped" builds cluster-respecting inner folds with lps.grouped.foldid() on the inner-training rows of cluster.id (passing the same per-fold offset seed when inner.shuffle.seed is supplied).

Value

A list of class "lps_nested_cv":

nested.rmse

pooled RMSE of the outer-test predictions over all rows with finite predictions.

n.missing.predictions

count of non-finite outer-test predictions (e.g. unstable.action = "na" fallbacks).

folds

one row per outer fold: fold label, test size, selected support.size / degree / kernel / bandwidth.multiplier, the inner selected-min CV score, and the fold's outer-test RMSE.

predictions

length-nrow(X) vector of outer-test predictions (each row predicted by the fold that held it out).

selected.min

the paired selected-min arm on the same outer.foldid: selected (the fit.lps selected row), cv.score (its observed CV RMSE), foldid, and fit (the full-data fit.lps object, usable for deployment on an external test set).

outer.foldid, test.index, train.index, inner.foldid, inner.foldid.used, inner.cv.table

complete fold/index telemetry: the realized inner fold ids as constructed and as recorded by each inner fit.lps ($foldid), per-fold index sets, and each fold's full inner CV candidate table.

outer.cluster.whole

TRUE/FALSE whether every cluster lies wholly inside one outer fold (NA when cluster.id is absent).

inner.folds, inner.foldid.method, inner.shuffle.seed, fit.args, call

the recorded configuration.

See Also

lps.grouped.foldid() for the grouped fold constructor.

Examples

X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
y <- sin(2 * pi * X[, 1])
nested <- lps.nested.cv(
  X, y, outer.foldid = rep(1:3, length.out = 18), inner.folds = 2L,
  fit.args = list(support.grid = 6L, degree.grid = 1L,
                  kernel.grid = "gaussian", backend = "R")
)
nested$rmse

Pointwise Variance and Confidence Band for a Fixed-Configuration LPS Fit

Description

Computes, from the analytically extracted linear-smoother matrix S of a fixed-configuration LPS fit (see smoother.matrix()), the pointwise variance ⁠Var(fitted_i) = sigma^2 * sum_j S_ij^2⁠, the effective degrees of freedom df = tr(S), the plug-in noise estimate sigma.hat^2 = RSS / (n - 2*tr(S) + sum(S^2)), and the pointwise confidence band ⁠fitted_i +/- z * sigma * ||S_i.||_2⁠ with z = qnorm(1 - (1 - level) / 2).

Usage

lps.pointwise.band(
  object,
  sigma = NULL,
  level = 0.95,
  check.tol = 1e-10,
  variance.method = c("residual", "legacy")
)

Arguments

Details

With sigma supplied (known noise standard deviation), the band uses it directly (sigma.source = "known"); with sigma = NULL the plug-in sigma.hat is used (sigma.source = "plug.in"). sigma.hat is reported in both modes. The variance is purely the linear-smoother sampling variance: the band makes no bias correction, so where bias dominates (boundary, high curvature) undercoverage is expected and must be reported, not masked.

Requires X.eval identical to X (the square training-fit smoother, on which tr(S) and RSS are defined), in addition to all restrictions of smoother.matrix(). If any evaluation point's local fit returned NA, its variance/band entries are NA, and df, rss, and sigma.hat are NA as well; supplying a known sigma still yields bands at the unaffected points.

Value

A list of class "lps.pointwise.band" with named fields: fitted (the fit's raw fitted values), se, variance, lower, upper, level, z, sigma (the supplied known sigma, or NA in plug-in mode), sigma.hat, sigma.source ("known" or "plug.in"), df (tr(S)), residual.df, variance.method, rss, n.train, smoother.row.norm (⁠||S_i.||_2⁠), and configuration (the pinned fit configuration).

See Also

smoother.matrix()

Examples

set.seed(1)
n <- 30
X <- matrix(runif(2 * n, -1, 1), ncol = 2)
y <- sin(pi * X[, 1]) + 0.1 * rnorm(n)
fit <- fit.lps(X, y, foldid = rep(1:2, length.out = n),
               support.grid = 12L, degree.grid = 1L,
               kernel.grid = "tricube", backend = "R",
               design.basis = "orthogonal.polynomial.drop",
               ridge.multiplier.grid = 0, ridge.condition.max = Inf,
               unstable.action = "na")
band.known <- lps.pointwise.band(fit, sigma = 0.1)
band.plugin <- lps.pointwise.band(fit)
band.known$df
band.plugin$sigma.hat

Compute MALPS Linear-Smoother Diagnostics

Description

Computes effective degrees of freedom, generalized cross-validation (GCV), and analytic leave-one-out residual diagnostics from a MALPS smoother matrix. These diagnostics are exact for fixed-support, fixed-weight linear MALPS fits. For cross-validated fits they are conditional on the selected support profile. Robust fits are rejected by default for the same reason described in smoother.matrix.

Usage

malps.gcv(
  object,
  y = NULL,
  smoother.matrix = NULL,
  include.loocv = TRUE,
  max.n = 1000L,
  allow.robust = FALSE,
  ...
)

Arguments

Details

The GCV score is computed as

\mathrm{GCV} = \frac{n^{-1}\|y - S y\|_2^2} {(1 - \operatorname{tr}(S)/n)^2}.

The analytic leave-one-out residuals are computed as

e_i^{\mathrm{loo}} = \frac{y_i - \hat y_i}{1 - S_{ii}}.

Value

A list with residual, fitted-value, EDF, GCV, and optional LOOCV diagnostics.

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
fit <- fit.malps(X, X[, 1]^2, degree = 1L,
                 support.type = "knn", support.size = 8L)
malps.gcv(fit)$gcv

Materialize a synthetic-data specification

Description

Materialize a synthetic-data specification

Usage

materialize.synthetic(
  spec,
  n = NULL,
  seed,
  rng.policy = c("named.stream.v1", "legacy"),
  validate = TRUE
)

Arguments

Value

A synthetic_dataset.

Examples

spec <- synthetic.registry.spec("G1")
materialize.synthetic(spec, n = 20L, seed = 1L)

Materialize a frozen synthetic instance

Description

Materialize a frozen synthetic instance

Usage

materialize.synthetic.instance(instance.id, validate = TRUE)

Arguments

Value

A synthetic_dataset whose dataset.id is instance.id.

Examples

materialize.synthetic.instance("geosmooth.g1.default.v1")

Construct a Graph Heat-Time Grid

Description

Creates a positive heat-time grid from a "metric.graph.lowpass.basis" object.

Usage

metric.graph.heat.eta.grid(
  basis,
  rule = c("spectral_guarded", "w1_inverse_spectrum"),
  n.initial = 40L,
  include.zero = FALSE,
  truncation.tol = 1e-04,
  equilibrium.tol = 1e-04
)

Arguments

Value

A numeric vector with grid-construction metadata stored as attributes. Exact null modes alone are excluded; uniformly rescaling conductances inversely rescales heat times. Unrepresentable endpoints cause an error.

Examples

adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
lengths <- list(1, c(1, 1), c(1, 1), 1)
basis <- metric.graph.lowpass.basis(adj, lengths, n.eigenpairs = 4L,
                                    eigen.solver = "dense")
metric.graph.heat.eta.grid(basis, n.initial = 6L)

Propose a Guarded Lower Graph Heat Time

Description

Proposes one lower positive heat time after an external search controller has classified the current lower endpoint. The helper enforces a geometric expansion step, an identity-departure floor, a maximum number of expansion rounds, and the truncated-basis resolution certificate. It does not decide whether an endpoint is competitive; fold-level or cohort-level search logic remains the caller's responsibility.

Usage

metric.graph.heat.extend.lower(
  basis,
  eta.grid,
  endpoint.status = c("active", "inactive", "unresolved"),
  expansion.factor = 3,
  max.expansions = 3L,
  expansions.completed = 0L,
  identity.departure = 0.01,
  truncation.tol = 1e-04,
  unresolved.action = c("error", "mark")
)

Arguments

Value

A list of class "metric.graph.heat.lower.extension" with the augmented or unchanged grid and proposal telemetry. A proposal is added only when admitted is TRUE.

Examples

adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
lengths <- list(1, c(1, 1), c(1, 1), 1)
basis <- metric.graph.lowpass.basis(adj, lengths, n.eigenpairs = 4L,
                                    eigen.solver = "dense")
eta <- metric.graph.heat.eta.grid(basis, n.initial = 6L)
metric.graph.heat.extend.lower(basis, eta, endpoint.status = "active")

Construct a Reusable Metric Graph Low-Pass Spectral Basis

Description

Builds a metric-conductance graph Laplacian and computes the low-frequency eigensystem used by graph low-pass filters. The resulting object is response-independent and can be reused for many responses and filter parameters.

Usage

metric.graph.lowpass.basis(
  adj.list,
  weight.list,
  conductance.rule = c("inverse.length.power", "exp.length", "exp.length.squared",
    "self.tuned.gaussian"),
  conductance.epsilon = 1e-08,
  conductance.alpha = 1,
  conductance.sigma = NULL,
  conductance.sigma.rule = c("edge.quantile", "median", "local.k"),
  conductance.sigma.quantile = 0.75,
  conductance.local.k = 5L,
  laplacian.type = c("unnormalized", "symmetric.normalized"),
  n.eigenpairs = 50L,
  eigen.solver = c("auto", "sparse", "dense"),
  dense.eigen.threshold = 200L,
  dense.fallback.threshold = 5000L,
  dense.fallback = c("auto", "never", "always"),
  verbose = FALSE,
  weight.type = c("length", "conductance"),
  spectrum = c("partial", "full"),
  numerical.policy = c("conservative", "verified.fallback")
)

Arguments

Details

Numerical verification anchors the known component nullspace (square-root degree for the normalized operator), and rejects unresolved positive modes rather than clamping them to zero. Heat paths, fits and refits additionally check spectral error at the requested times. Missing-response blocks must pass a conditioning check as well as Cholesky and residual checks; failure returns an error without deleting positive edges or adding a ridge. These are conservative double-precision diagnostics, not rigorous interval error bounds. A complete basis alone does not certify arbitrary long times.

The basis is complete only when it contains one eigenvector per graph vertex. A truncated basis represents only the retained low-frequency subspace. The largest retained eigenvalue supplies a conservative proxy for bounding the contribution of omitted modes because all omitted eigenvalues are at least as large.

Value

A list of class "metric.graph.lowpass.basis" containing the graph operator, eigenvalues, eigenvectors, solver metadata, and spectral completeness diagnostics.

Examples

adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
lengths <- list(1, c(1, 1), c(1, 1), 1)
metric.graph.lowpass.basis(adj, lengths, n.eigenpairs = 4L,
                           eigen.solver = "dense")

Construct a Metric-Conductance Graph Low-Pass Operator

Description

Builds a weighted graph Laplacian from supplied conductances or by transforming metric edge lengths into conductances. This operator is a direct metric-conductance comparator to the legacy rdgraph regression smoother, not a replacement for the Riemannian-complex/overlap-density smoother.

Usage

metric.graph.lowpass.operator(
  adj.list,
  weight.list,
  conductance.rule = c("inverse.length.power", "exp.length", "exp.length.squared",
    "self.tuned.gaussian"),
  conductance.epsilon = 1e-08,
  conductance.alpha = 1,
  conductance.sigma = NULL,
  conductance.sigma.rule = c("edge.quantile", "median", "local.k"),
  conductance.sigma.quantile = 0.75,
  conductance.local.k = 5L,
  laplacian.type = c("unnormalized", "symmetric.normalized"),
  return.sparse = TRUE,
  verbose = FALSE,
  weight.type = c("length", "conductance")
)

Arguments

Details

The legacy rdgraph regression precomputed-graph path uses supplied weight.list values as edge lengths for neighborhood ordering. The spectral conductance in that smoother is overlap-density based:

c_e^\rho = 1 / \max(\rho_1(e), 10^{-10}),

where \rho_1(e) is computed by the Riemannian-complex density machinery.

This function instead constructs conductances directly from metric lengths:

c_{ij} = \phi(\ell_{ij}).

Supported transforms are

c_{ij}=(\ell_{ij}+\epsilon)^{-\alpha},

c_{ij}=\exp(-\ell_{ij}/\sigma),

c_{ij}=\exp(-\ell_{ij}^{2}/\sigma^{2}),

and

c_{ij}=\exp(-\ell_{ij}^{2}/(\sigma_i\sigma_j)).

For laplacian.type = "symmetric.normalized", smoothing is performed in the Euclidean eigenbasis of L_{\mathrm{sym}}. The null vector is proportional to \sqrt{d}, not the constant vector, so this mode does not preserve constant responses in the same way as the unnormalized Laplacian.

Value

A list of class "metric.graph.lowpass.operator" containing the edge table, conductances, degree vector, Laplacian triplets, summaries, and optionally a sparse Laplacian matrix.

Examples

adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
lengths <- list(1, c(1, 1), c(1, 1), 1)
metric.graph.lowpass.operator(adj, lengths)

Summarize a Geometric Smoother or Probability-Mass Fit

Description

Summaries expose stored settings and diagnostics without recomputing geometry, fitting, or selection. Console output is bounded; the summary's structured fields retain selected controls, criterion tables, local diagnostic counts, solver information, and density accounting for closer inspection.

Usage

## S3 method for class 'summary.geosmooth_fit'
print(x, ...)

## S3 method for class 'lps'
summary(object, ...)

## S3 method for class 'malps'
summary(object, ...)

## S3 method for class 'ps_lps'
summary(object, ...)

## S3 method for class 'chart_kernel'
summary(object, ...)

## S3 method for class 'local_likelihood'
summary(object, ...)

## S3 method for class 'lpl_tf'
summary(object, ...)

## S3 method for class 'slpl_tf'
summary(object, ...)

## S3 method for class 'metric.graph.lowpass.fit'
summary(object, ...)

## S3 method for class 'metric.graph.lowpass.refit'
summary(object, ...)

## S3 method for class 'ssrhe.hessian.fit'
summary(object, ...)

## S3 method for class 'ssrhe.hessian.refit'
summary(object, ...)

## S3 method for class 'ssrhe.hessian.l1.fit'
summary(object, ...)

## S3 method for class 'graph.trend.filtering.fit'
summary(object, ...)

## S3 method for class 'pttf.trend.filtering.fit'
summary(object, ...)

## S3 method for class 'density_fit'
summary(object, ...)

## S3 method for class 'lpl_tf'
print(x, ...)

## S3 method for class 'slpl_tf'
print(x, ...)

## S3 method for class 'density_fit'
print(x, ...)

## S3 method for class 'ps_lps'
print(x, ...)

## S3 method for class 'chart_kernel'
print(x, ...)

## S3 method for class 'local_likelihood'
print(x, ...)

Arguments

Details

A finite result is not a convergence certificate. Where the algorithm records no stopping certificate, convergence is explicitly unavailable. Local fallback counts and failed evaluations are reported separately. Selection criteria are conditional tuning scores; they are not an independent assessment of prediction error. GCV, CV MSE, RMSE, and likelihood losses have different interpretations.

Value

summary() returns a summary.geosmooth_fit list with model, n.training, n.evaluation, n.responses, n.observed, nonfinite, status, convergence, prediction.domain, selected, criteria, diagnostics, solver, mass, accounting, and warnings fields. Print methods return their input invisibly.

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
fit <- fit.chart.kernel(X, X[, 1]^2, support.size = 8L)
fit
summary(fit)
mass <- normalize.density(c(1, 2, 3))
mass
summary(mass)$accounting

Extract Stored Fitted Values and Response Residuals

Description

Use fitted(object) for the stored response-scale fitted values, preserving vector or matrix shape, row order, names, and missing values. It does not recompute a fit or imply prediction at new points. For harmonic.smoother() results it returns harmonic_predictions. The basic perform.harmonic.smoothing() result remains an unclassed list, accessed through ⁠$harmonic_predictions⁠; density probability masses remain ⁠$rho⁠.

Usage

## S3 method for class 'lps'
fitted(object, ...)

## S3 method for class 'lps'
residuals(object, ...)

## S3 method for class 'malps'
fitted(object, ...)

## S3 method for class 'malps'
residuals(object, ...)

## S3 method for class 'ps_lps'
fitted(object, ...)

## S3 method for class 'ps_lps'
residuals(object, ...)

## S3 method for class 'chart_kernel'
fitted(object, ...)

## S3 method for class 'chart_kernel'
residuals(object, ...)

## S3 method for class 'local_likelihood'
fitted(object, ...)

## S3 method for class 'local_likelihood'
residuals(object, ...)

## S3 method for class 'lpl_tf'
fitted(object, ...)

## S3 method for class 'lpl_tf'
residuals(object, ...)

## S3 method for class 'slpl_tf'
fitted(object, ...)

## S3 method for class 'slpl_tf'
residuals(object, ...)

## S3 method for class 'metric.graph.lowpass.fit'
fitted(object, ...)

## S3 method for class 'metric.graph.lowpass.fit'
residuals(object, ...)

## S3 method for class 'metric.graph.lowpass.refit'
fitted(object, ...)

## S3 method for class 'metric.graph.lowpass.refit'
residuals(object, ...)

## S3 method for class 'ssrhe.hessian.fit'
fitted(object, ...)

## S3 method for class 'ssrhe.hessian.fit'
residuals(object, ...)

## S3 method for class 'ssrhe.hessian.refit'
fitted(object, ...)

## S3 method for class 'ssrhe.hessian.refit'
residuals(object, ...)

## S3 method for class 'ssrhe.hessian.l1.fit'
fitted(object, ...)

## S3 method for class 'ssrhe.hessian.l1.fit'
residuals(object, ...)

## S3 method for class 'graph.trend.filtering.fit'
fitted(object, ...)

## S3 method for class 'graph.trend.filtering.fit'
residuals(object, ...)

## S3 method for class 'pttf.trend.filtering.fit'
fitted(object, ...)

## S3 method for class 'pttf.trend.filtering.fit'
residuals(object, ...)

## S3 method for class 'harmonic_smoother'
fitted(object, ...)

Arguments

Details

residuals(object) returns stored response residuals where available, including missing-label masks. Otherwise it computes y - fitted.values only when the evaluation coordinates equal the training coordinates. Residuals at different evaluation points are undefined and cause an error, even when training and evaluation have the same number of rows. Binomial residuals are on the response (probability) scale, not Pearson or deviance residuals. These accessors do not change which refit results are reusable.

Value

The stored vector or matrix, or response residuals of the same shape.

See Also

refit, smoother.matrix

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
model <- fit.malps(X, X[, 1]^2, degree = 1L,
                   support.type = "knn", support.size = 8L)
head(fitted(model))
head(residuals(model))

Normalize Fitted Values Into a Density

Description

Converts a numeric field or fitted smoother/regression object into a probability mass vector over the evaluation support. This is the explicit adapter between ordinary smoothers such as fit.lps, fit.ps.lps, and fit.metric.graph.lowpass and density or occupation-density workflows.

Usage

normalize.density(x, ...)

## S3 method for class 'numeric'
normalize.density(
  x,
  X = NULL,
  density.control = list(),
  method.id = "normalized_numeric",
  keep.source.fit = FALSE,
  adj.list = NULL,
  empirical.rho = NULL,
  return.details = TRUE,
  ...
)

## Default S3 method:
normalize.density(
  x,
  X = NULL,
  density.control = list(),
  method.id = NULL,
  keep.source.fit = TRUE,
  adj.list = NULL,
  empirical.rho = NULL,
  return.details = TRUE,
  ...
)

## S3 method for class 'lps'
normalize.density(
  x,
  X = NULL,
  density.control = list(),
  method.id = NULL,
  keep.source.fit = TRUE,
  adj.list = NULL,
  empirical.rho = NULL,
  return.details = TRUE,
  ...
)

## S3 method for class 'ps_lps'
normalize.density(
  x,
  X = NULL,
  density.control = list(),
  method.id = NULL,
  keep.source.fit = TRUE,
  adj.list = NULL,
  empirical.rho = NULL,
  return.details = TRUE,
  ...
)

## S3 method for class 'metric.graph.lowpass.fit'
normalize.density(
  x,
  X = NULL,
  density.control = list(),
  method.id = NULL,
  keep.source.fit = TRUE,
  adj.list = NULL,
  empirical.rho = NULL,
  return.details = TRUE,
  ...
)

## S3 method for class 'metric.graph.lowpass.refit'
normalize.density(
  x,
  X = NULL,
  density.control = list(),
  method.id = NULL,
  keep.source.fit = TRUE,
  adj.list = NULL,
  empirical.rho = NULL,
  return.details = TRUE,
  ...
)

Arguments

Value

A list of class "density_fit".

Examples

normalize.density(c(0.5, -0.1, 0.6, 0), X = matrix(1:4, ncol = 1))

Perform Harmonic Smoothing on Graph Function Values

Description

Applies harmonic smoothing to function values defined on vertices of a graph, preserving values at the boundary of a specified region while smoothly interpolating interior values. This function implements a discrete Laplace equation solution using weighted averaging.

Usage

perform.harmonic.smoothing(
  adj.list,
  weight.list,
  values,
  region.vertices,
  max.iterations = 100,
  tolerance = 1e-06
)

Arguments

Details

Boundary vertices are exactly the region vertices with a neighbor outside the region. Degree-one vertices are not automatically boundary vertices. Boundary values and values outside the region remain fixed. Each synchronous relaxation step replaces interior values by neighbor averages with conductance 1 / (w + 10^{-10}), where w is the supplied edge length.

With no boundary, this function returns the input unchanged, with status "no_boundary" and converged = FALSE. With no interior there is nothing to solve: status is "no_interior" and convergence is true. On a disconnected region, components without a connection to fixed boundary values are not uniquely anchored; relaxation can oscillate. Inspect convergence before using the result as a harmonic solution. Unlike this basic function, harmonic.smoother() also runs relaxation when the whole region has no boundary, retaining its historical behavior.

Value

A list containing:

• harmonic_predictions: the numeric vector of smoothed values;

• converged: whether both the maximum update and discrete harmonic residual fell below tolerance;

• num_region, num_boundary, and num_interior: vertex counts for the requested region and its boundary/interior split;

• num_iterations: the number of relaxation iterations run;

• max_change and max_residual: per-iteration convergence diagnostics;

• status: converged, iteration_limit, no_boundary, or no_interior.

See Also

harmonic.smoother for smoothing with topology tracking, get.region.boundary for boundary vertex identification

Examples

adj.list <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
weight.list <- list(1, c(1, 1), c(1, 1), 1)
values <- c(0, 2, -1, 1)
result <- perform.harmonic.smoothing(
  adj.list, weight.list, values, region.vertices = 1:3
)
result$harmonic_predictions

Plot Method for Harmonic Smoother Results

Description

Creates various plots to visualize the results of harmonic smoothing with topology tracking.

Usage

## S3 method for class 'harmonic_smoother'
plot(x, y = NULL, ..., type = c("topology", "extrema", "values"))

Arguments

Value

Invisibly returns the input object.

See Also

harmonic.smoother

Examples

adj.list <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
weight.list <- list(1, c(1, 1), c(1, 1), 1)
values <- c(0, 2, -1, 1)
result <- harmonic.smoother(
  adj.list, weight.list, values, region.vertices = seq_along(values),
  max.iterations = 20, record.frequency = 2
)
plot(result, type = "values")

Plot a canonical synthetic dataset

Description

The default view = "auto" preserves the dimension-dependent display: one latent coordinate shows response and truth against that coordinate; otherwise the first two predictor coordinates show a geometry projection. Select "geometry" explicitly to see a curve's embedding, or "response" to see responses against the first latent coordinate (observation index when latent coordinates are absent). A two-coordinate projection can hide separation in higher ambient dimensions. Plotting does not draw random values.

Usage

## S3 method for class 'synthetic_dataset'
plot(
  x,
  color = c("response", "truth", "region"),
  ...,
  view = c("auto", "geometry", "response"),
  legend = TRUE,
  col = NULL
)

Arguments

Value

x, invisibly.

Examples

x <- materialize.synthetic(synthetic.registry.spec("G1"), n = 20L, seed = 1L)
plot(x)
plot(x, view = "geometry", color = "truth")

Predict from Gaussian Trend Filtering

Description

Predict from Gaussian Trend Filtering

Usage

## S3 method for class 'gaussian.trend.filtering.fit'
predict(
  object,
  newdata = NULL,
  covariates = NULL,
  offset = NULL,
  type = c("response", "link"),
  lambda.index = NULL,
  allow.unconverged = FALSE,
  ...
)

## S3 method for class 'gaussian.trend.filtering.fit'
fitted(object, ...)

## S3 method for class 'gaussian.trend.filtering.fit'
residuals(object, ...)

## S3 method for class 'gaussian.trend.filtering.fit'
print(x, ...)

Arguments

Value

A numeric vector for one path index, otherwise a matrix. Attributes outside.range and lambda record boundary predictions and penalties.

Examples

x <- seq(0, 1, length.out = 20)
fit <- fit.gaussian.trend.filtering(x, 1+x-x^2, lambda = Inf)
predict(fit, c(0.25, 0.75))
head(fitted(fit))
head(residuals(fit))
print(fit)

Predict from Logistic Trend Filtering

Description

Predict from Logistic Trend Filtering

Usage

## S3 method for class 'logistic.trend.filtering.fit'
predict(
  object,
  newdata = NULL,
  covariates = NULL,
  offset = NULL,
  type = c("response", "link"),
  lambda.index = NULL,
  allow.unconverged = FALSE,
  ...
)

## S3 method for class 'logistic.trend.filtering.fit'
fitted(object, ...)

## S3 method for class 'logistic.trend.filtering.fit'
residuals(object, ...)

## S3 method for class 'logistic.trend.filtering.fit'
print(x, ...)

Arguments

Value

A numeric vector for one path index, otherwise a matrix. Attributes outside.range and lambda identify boundary predictions and penalties.

Examples

x <- seq(0, 1, length.out = 24)
fit <- fit.logistic.trend.filtering(x, rep(c(0, 1), 12), lambda = Inf)
predict(fit, c(0.25, 0.75))
head(fitted(fit))
head(residuals(fit))
print(fit)

Predict From A Local Polynomial Lifting Trend Filter

Description

Predict From A Local Polynomial Lifting Trend Filter

Usage

## S3 method for class 'lpl_tf'
predict(
  object,
  newdata = NULL,
  type = c("response"),
  allow.incomplete = FALSE,
  ...
)

Arguments

Value

Fitted values at the training points.

Examples

if (requireNamespace("genlasso", quietly = TRUE)) {
  X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
  fit <- fit.lpl.tf(X, sin(2 * pi * X[, 1]), degree = 1L,
                    support.type = "knn", support.size = 7L,
                    lambda = 0.1, lambda.selection = "fixed")
  predict(fit)
}

Predict from an LPS Fit

Description

Predict from an LPS Fit

Usage

## S3 method for class 'lps'
predict(object, newdata = NULL, type = c("response", "raw"), ...)

Arguments

Value

A numeric vector of predictions.

Examples

X <- matrix(seq(0, 1, length.out = 12), ncol = 1)
fit <- fit.lps(X, sin(2 * pi * X[, 1]), support.grid = 6L)
predict(fit)

Predict From A MALPS Fit

Description

Predict From A MALPS Fit

Usage

## S3 method for class 'malps'
predict(
  object,
  newdata = NULL,
  type = c("response"),
  allow.incomplete = FALSE,
  ...
)

Arguments

Value

Numeric vector of fitted values.

Examples

X <- matrix(seq(0, 1, length.out = 12), ncol = 1)
fit <- fit.malps(X, sin(2 * pi * X[, 1]), support.size = 6L)
predict(fit)

Predict Nearest-Person Binary Probabilities

Description

Predict Nearest-Person Binary Probabilities

Usage

## S3 method for class 'neighbor.probability.fit'
predict(
  object,
  newdata = NULL,
  distances = NULL,
  query.group = NULL,
  k = object$k,
  type = c("response", "diagnostics"),
  ...
)

Arguments

Value

A query-by-size matrix, or a diagnostic list with that matrix as probability. Unsupported estimates remain NA. boundary.ties counts distances within 1e-12 of the cutoff; ordering itself uses exact distance ties. Diagnoses count distinct people. With positive eps, boundary counts cover retrieved participants only (boundary.complete = FALSE). search.candidates gives the final number of candidate records, including visit refinement.

Examples

X <- matrix(0:5, ncol = 1)
f <- fit.neighbor.probability(X, c(1, 0, 1, 0, 0, 1), k = 2)
predict(f) # leave-participant-out, not resubstitution
predict(f, matrix(2.5, ncol = 1), type = "diagnostics")

Predict From A Fixed-Operator Synchronized LPL-TF Model

Description

Predict From A Fixed-Operator Synchronized LPL-TF Model

Usage

## S3 method for class 'slpl_tf'
predict(object, newdata = NULL, type = c("response"), ...)

Arguments

Value

Fitted values at the training points.

Examples

if (requireNamespace("genlasso", quietly = TRUE)) {
  X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
  fit <- fit.slpl.tf(X, sin(2 * pi * X[, 1]), degree = 1L,
                     support.type = "knn", support.size = 7L,
                     lambda1 = 0.1, lambda2 = 0,
                     lambda.selection = "fixed")
  predict(fit)
}

Print Geosmooth S3 Objects

Description

Print concise summaries for fitted geosmooth models and constructed operators.

Usage

## S3 method for class 'graph.trend.filtering.fit'
print(x, ...)

## S3 method for class 'lpl_tf_operator'
print(x, ...)

## S3 method for class 'lps'
print(x, ...)

## S3 method for class 'malps'
print(x, ...)

## S3 method for class 'malps_bootstrap'
print(x, ...)

## S3 method for class 'metric.graph.lowpass.fit'
print(x, ...)

## S3 method for class 'metric.graph.lowpass.refit'
print(x, ...)

## S3 method for class 'pttf.trend.filtering.fit'
print(x, ...)

## S3 method for class 'pttf_operator'
print(x, ...)

## S3 method for class 'slpl_tf_operator'
print(x, ...)

## S3 method for class 'ssrhe.hessian.operator'
print(x, ...)

## S3 method for class 'ssrhe.hessian.l1.fit'
print(x, ...)

## S3 method for class 'ssrhe.hessian.fit'
print(x, ...)

## S3 method for class 'ssrhe.hessian.refit'
print(x, ...)

## S3 method for class 'ssrhe.hessian.cv.fit'
print(x, ...)

## S3 method for class 'ssrhe.hessian.gcv.fit'
print(x, ...)

## S3 method for class 'transported.graph.hessian.operator'
print(x, ...)

Arguments

Value

Invisibly returns x.

Examples

X <- matrix(seq(0, 1, length.out = 12), ncol = 1)
fit <- fit.lps(X, sin(2 * pi * X[, 1]), support.grid = 6L)
print(fit)

Print Method for Harmonic Smoother Results

Description

Prints a concise summary of harmonic smoother results.

Usage

## S3 method for class 'harmonic_smoother'
print(x, ...)

Arguments

Value

Invisibly returns the input object.

See Also

harmonic.smoother, summary.harmonic_smoother

Examples

fit <- structure(
  list(
    stable_iteration = 2L,
    i_harmonic_predictions = matrix(1:4, nrow = 2),
    i_basins = list(matrix(c(1, 1), ncol = 2)),
    topology_differences = numeric()
  ),
  class = "harmonic_smoother"
)
print(fit)

Print Method for Summary of Harmonic Smoother Results

Description

Prints the summary of harmonic smoother results in a formatted manner.

Usage

## S3 method for class 'summary.harmonic_smoother'
print(x, ...)

Arguments

Value

Invisibly returns the input object.

See Also

summary.harmonic_smoother

Examples

object <- structure(
  list(
    stable_iteration = 2L,
    iterations_recorded = 1L,
    initial_extrema = 1L,
    initial_maxima = 1L,
    initial_minima = 0L,
    final_extrema = 1L,
    final_maxima = 1L,
    final_minima = 0L,
    extrema_reduction = 0L,
    topology_diff_summary = NULL
  ),
  class = "summary.harmonic_smoother"
)
print(object)

Print a canonical synthetic dataset

Description

Print a canonical synthetic dataset

Usage

## S3 method for class 'synthetic_dataset'
print(x, ...)

Arguments

Value

x, invisibly.

Examples

x <- materialize.synthetic(synthetic.registry.spec("G1"), n = 20L, seed = 1L)
print(x)

Construct PTTF Local Geometry

Description

Builds the geometry layer for the experimental parallel-transport trend filtering (PTTF) route. The returned object contains graph supports, local PCA tangent frames, edge directions, orthogonal edge transports, sampling density metadata, and diagnostics. It does not assemble transported difference operators and it does not fit a response.

Usage

pttf.geometry(
  X,
  adj.list = NULL,
  weight.list = NULL,
  graph = c("rknn", "supplied"),
  tangent.dim,
  local.support = c("graph.disk", "neighbors"),
  min.support = NULL,
  max.hops = 3L,
  support.k = NULL,
  graph.k.scale = 1L,
  graph.radius.factor = 1,
  graph.radius.rule = "geomean",
  transport.rule = c("procrustes"),
  synchronize.orientation = TRUE,
  density.method = c("support.radius"),
  alpha = 1,
  diagnostics = TRUE
)

Arguments

Details

The stored transport convention is column-vector based. For an oriented edge i \leftarrow j, Oij maps tangent coordinates from frame j into frame i:

v_j \mapsto O_{ij} v_j.

For shared-support Procrustes, row-coordinate matrices Z_i and Z_j are formed on the same vertex set. The row problem \min_Q \|Z_j Q - Z_i\|_F^2 gives Q = U V^\top when Z_j^\top Z_i = U\Sigma V^\top. The stored column map is O_{ij}=Q^\top=VU^\top; the reported residual is \|Z_j O_{ij}^\top - Z_i\|_F / \max(\|Z_i\|_F,\epsilon).

Value

A list of class "pttf_geometry".

Examples

n <- 10L
X <- matrix(seq(0, 1, length.out = n), ncol = 1)
adj <- lapply(seq_len(n), function(i) intersect(c(i - 1L, i + 1L), 1:n))
lengths <- Map(function(i, j) abs(X[j, 1] - X[i, 1]), seq_len(n), adj)
pttf.geometry(X, adj, lengths, graph = "supplied", tangent.dim = 1L)

Construct PTTF Transported-Difference Operators

Description

Assembles experimental parallel-transport trend-filtering (PTTF) transported-difference operators from a pttf.geometry object. This is an operator-construction function only: it does not fit a response, tune a penalty parameter, or run an \ell_1/\ell_2 smoother.

Usage

pttf.operator(
  geometry,
  derivative.order = 3L,
  edge.status.policy = c("ok.only", "frame.fallback"),
  regression.weight.rule = c("inverse.length.squared", "inverse.length", "unit"),
  edge.length.epsilon = 1e-08,
  tensor.scaling = c("hs", "raw"),
  row.mass.rule = c("node.mass", "none"),
  row.normalize = c("none", "l2"),
  min.operator.rank.tol = 1e-10,
  max.operator.condition = 1e+08,
  return.intermediate = TRUE,
  return.A = TRUE,
  return.B = TRUE,
  diagnostics = TRUE
)

Arguments

Details

The operator keeps a full logical vertex/component layout internally. For order r, the logical field has one block of q_r=\binom{m+r-1}{r} components at each vertex. The returned A is compact: it includes rows only for accepted final-order vertex/component blocks, and row.table$full.row records the original full logical row.

For derivative level r, the local regression uses blocks \Phi_r(u_{ij})\in\mathbb R^{q_{r-1}\times q_r}, lifted weights W_i^{(r)}=\operatorname{diag}(w_{ij})\otimes I_{q_{r-1}}, and transported predecessor differences K_{r-1}(O_{ij})z_j-z_i. Unavailable predecessor blocks are never filled with zero; affected neighbors are dropped and recorded.

Value

A list of class "pttf_operator" containing sparse operators, row/vertex provenance, tensor-basis metadata, and diagnostics.

Examples

n <- 10L
X <- matrix(seq(0, 1, length.out = n), ncol = 1)
adj <- lapply(seq_len(n), function(i) intersect(c(i - 1L, i + 1L), 1:n))
lengths <- Map(function(i, j) abs(X[j, 1] - X[i, 1]), seq_len(n), adj)
geometry <- pttf.geometry(
  X, adj, lengths, graph = "supplied", tangent.dim = 1L
)
pttf.operator(geometry, derivative.order = 2L)

Filter Rows Of A PTTF Operator

Description

Creates a new "pttf_operator" object with a subset of compact operator rows. The sparse triplet payload is rebuilt from the filtered matrix so triplet row indices always refer to the filtered fit rows, while original compact and full logical row identities remain in row.table.

Usage

pttf.operator.filter.rows(
  operator,
  rows,
  reason = NULL,
  preserve.original = TRUE
)

Arguments

Value

A filtered "pttf_operator" object.

Examples

n <- 10L
X <- matrix(seq(0, 1, length.out = n), ncol = 1)
adj <- lapply(seq_len(n), function(i) intersect(c(i - 1L, i + 1L), 1:n))
lengths <- Map(function(i, j) abs(X[j, 1] - X[i, 1]), seq_len(n), adj)
geometry <- pttf.geometry(
  X, adj, lengths, graph = "supplied", tangent.dim = 1L
)
operator <- pttf.operator(geometry, derivative.order = 2L)
pttf.operator.filter.rows(operator, seq_len(max(1L, nrow(operator$A) - 2L)))

Refit A Model Using Its Stored Geometry

Description

Fit a new response while reusing the geometry or operators in a supported fitted object. R selects the method from the class of object.

Usage

refit(object, y, ...)

## Default S3 method:
refit(object, y, ...)

Arguments

Details

Refitting does not automatically repeat the original model selection:

malps

refit.malps reuses supports and averaging weights; accepts observation weights. Robust local residual weights, if used, are recomputed for the new response.

lpl_tf

refit.lpl_tf reuses the residual operator and penalty, or accepts a new fixed lambda.

slpl_tf

refit.slpl_tf reuses both operators and penalties, or accepts fixed lambda1 and lambda2.

metric.graph.lowpass.fit

refit.metric.graph.lowpass.fit reuses the spectral basis and smoothing parameter. Set per.column.gcv = TRUE to select a parameter for each response column.

ssrhe.hessian.fit

refit.ssrhe.hessian.fit reuses the Hessian operator and fixed penalties, with optional lambda1, lambda2, weights, and ridge. This also handles CV- and GCV-selected fits, retaining their selected penalties by default.

ssrhe.hessian.l1.fit

refit.ssrhe.hessian.l1.fit reuses the derivative operator with fixed or CV-selected penalties; accepts lambda.grid, lambda.selection, and solver controls.

Use the original fitted object for repeated metric graph low-pass or quadratic Hessian refits: their refit results are summaries, not reusable input fits. Other supported families return reusable fitted objects.

Value

A refitted object or refit summary, as documented by the method.

Migration

The six former exported refit.*() functions are replaced by refit(object, y, ...). Replace named fitted.model and y.new arguments with object and y, respectively. Method-specific controls and return classes are preserved. Call the generic rather than an individual S3 method; the former help names remain available for navigation.

See Also

smoother.matrix

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
fit <- fit.malps(X, X[, 1]^2, degree = 1L,
                 support.type = "knn", support.size = 8L)
updated <- refit(fit, y = X[, 1]^3)
head(predict(updated))

Refit A Local Polynomial Lifting Trend Filter

Description

Reuses the fixed operator and selected lambda from an existing LPL-TF fit.

Usage

## S3 method for class 'lpl_tf'
refit(object, y, lambda = NULL, reuse.lambda = TRUE, verbose = FALSE, ...)

Arguments

Value

A refitted "lpl_tf" object.

See Also

refit for supported fitted objects and migration.

Examples

if (requireNamespace("genlasso", quietly = TRUE)) {
  X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
  fit <- fit.lpl.tf(X, X[, 1]^2, degree = 1L,
                    support.type = "knn", support.size = 7L,
                    lambda = 0.1, lambda.selection = "fixed")
  refit(fit, y = X[, 1]^3)
}

Refit A MALPS Smoother With A New Response

Description

Reuses the anchors, selected support profile, local supports, prediction supports, coordinate method, kernel, degree, and local-solver controls from a previous fit.malps result, then recomputes local polynomial coefficients for a new response vector. Observation weights, when supplied, multiply the stored geometric fitting weights for local coefficient estimation. Supports and model-averaging weights remain fixed. A weighted refit requires every stored local support to retain at least one positive effective fitting weight; sparse bootstrap-style weights can therefore fail even when the global weight vector is not all zero. If the original fit used robust local residual reweighting, robust weights are recomputed for the new response because they are response-dependent.

Usage

## S3 method for class 'malps'
refit(
  object,
  y = NULL,
  weights = NULL,
  reuse.selection = TRUE,
  refit.local.coefficients = TRUE,
  verbose = FALSE,
  ...
)

Arguments

Value

A "malps" object with the same support profile and new fitted values.

See Also

refit for supported fitted objects and migration.

Examples

X <- matrix(seq(0, 1, length.out = 20), ncol = 1)
fit <- fit.malps(X, X[, 1]^2, degree = 1L,
                 support.type = "knn", support.size = 8L)
refit(fit, y = X[, 1]^3)

Refit Metric-Conductance Graph Low-Pass Regression

Description

Reuses a fitted metric graph low-pass eigensystem to smooth new responses.

Usage

## S3 method for class 'metric.graph.lowpass.fit'
refit(
  object,
  y,
  per.column.gcv = FALSE,
  eta.grid = NULL,
  n.candidates = 40L,
  n.cores = 1L,
  block.size = NULL,
  verbose = FALSE,
  numerical.policy = NULL,
  ...
)

Arguments

Value

A list of class "metric.graph.lowpass.refit".

See Also

refit for supported fitted objects and migration.

Examples

adj <- list(2L, c(1L, 3L), c(2L, 4L), 3L)
lengths <- list(1, c(1, 1), c(1, 1), 1)
fit <- fit.metric.graph.lowpass(
  adj, lengths, y = 1:4, n.eigenpairs = 4L,
  eta.grid = c(0.1, 1), eigen.solver = "dense"
)
refit(fit, y = 4:1)

Refit A Fixed-Operator Synchronized LPL-TF Model

Description

Refit A Fixed-Operator Synchronized LPL-TF Model

Usage

## S3 method for class 'slpl_tf'
refit(
  object,
  y,
  lambda1 = NULL,
  lambda2 = NULL,
  reuse.lambda = TRUE,
  verbose = FALSE,
  ...
)

Arguments

Value

A refitted "slpl_tf" object.

See Also

refit for supported fitted objects and migration.

Examples

if (requireNamespace("genlasso", quietly = TRUE)) {
  X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
  fit <- fit.slpl.tf(X, X[, 1]^2, degree = 1L,
                     support.type = "knn", support.size = 7L,
                     lambda1 = 0.1, lambda2 = 0,
                     lambda.selection = "fixed")
  refit(fit, y = X[, 1]^3)
}

Refit SSRHE-Style Hessian-Energy Regression

Description

Reuses the operator from fit.ssrhe.hessian.regression to fit new responses or new fixed penalty weights without rebuilding local PCA neighborhoods or Hessian-energy matrices.

Usage

## S3 method for class 'ssrhe.hessian.fit'
refit(
  object,
  y = NULL,
  lambda1 = object$lambda$lambda1,
  lambda2 = object$lambda$lambda2,
  weights = NULL,
  ridge = object$lambda$ridge,
  verbose = FALSE,
  ...
)

Arguments

Value

A list of class "ssrhe.hessian.refit".

See Also

refit for supported fitted objects and migration.

Examples

X <- matrix(seq(0, 1, length.out = 12), ncol = 1)
fit <- fit.ssrhe.hessian.regression(
  X, X[, 1]^2, k = 6L, tangent.dim = 1L, lambda1 = 0.1
)
refit(fit, y = X[, 1]^3)

Refit SSRHE-Style Hessian L1 Regression

Description

Reuses the A operator from fit.ssrhe.hessian.l1.regression to fit a new response or a new lambda grid without rebuilding local neighborhoods or local PCA Hessian rows.

Usage

## S3 method for class 'ssrhe.hessian.l1.fit'
refit(
  object,
  y = NULL,
  lambda.grid = object$lambda,
  lambda.selection = c("fixed", "cv"),
  weights = NULL,
  n.lambda = 40L,
  nfolds = 5L,
  fold.id = NULL,
  loss = c("mse", "mae"),
  selection = c("min", "one.se"),
  solver = c("genlasso", "admm", "auto"),
  row.scaling = c("none", "l2"),
  admm.rho = NULL,
  admm.maxiter = 2000L,
  admm.abstol = 1e-06,
  admm.reltol = 1e-05,
  maxsteps = 2000L,
  minlam = 0,
  approx = FALSE,
  rtol = 1e-07,
  btol = 1e-07,
  eps = 1e-04,
  verbose = FALSE,
  admm.adaptive.rho = TRUE,
  ...
)

Arguments

Value

A list of class "ssrhe.hessian.l1.refit".

See Also

refit for supported fitted objects and migration.

Examples

X <- matrix(seq(0, 1, length.out = 12), ncol = 1)
fit <- fit.ssrhe.hessian.l1.regression(
  X, X[, 1]^2, k = 6L, tangent.dim = 1L,
  lambda.grid = 0.05, lambda.selection = "fixed", solver = "admm"
)
refit(
  fit, y = X[, 1]^3, solver = "admm"
)

Evaluate Nested Selection Among Candidate Graphs

Description

Chooses the graph using inner-fold loss only, then assembles its outer-fold predictions. Only graphs with full inner-fold prediction coverage compete; unsupported outer predictions remain NA. Candidate vertex identities and ordering must agree. Final graph selection uses full-data inner CV.

Usage

select.masked.graph.field(comparison, study)

Arguments

Value

Selected graph per outer fold, outer predictions and coverage, participant-balanced MSE, final graph name, and its final fitted field.

Examples

op <- metric.graph.lowpass.operator(list(2L, c(1L, 3L), 2L),
                                   list(1, c(1, 1), 1))
studies <- list(example = list(vertex = rep(c(1, 3), 6), y = 1:12,
                               group = letters[1:12]))
comparison <- compare.masked.graph.fields(list(path = op), studies,
    lambdas = c(0.01, 1), outer.folds = 3, inner.folds = 2)
selected <- select.masked.graph.field(comparison, "example")
selected$final.graph
selected$coverage

Construct A Synchronized LPL-TF Fixed Operator

Description

Builds the fixed S-LPL-TF operator pair. The LPL component A_{\mathrm{LPL}} is the accepted self-excluded LPL-TF residual operator from lpl.tf.operator. The synchronization component C_{\mathrm{sync}} compares inclusive local polynomial prediction maps over overlapping supports.

Usage

slpl.tf.operator(
  X,
  adj.list = NULL,
  weight.list = NULL,
  graph = NULL,
  graph.stage = "final",
  anchor.index = NULL,
  anchor.coordinates = NULL,
  degree = 2L,
  support.type = c("adaptive.radius", "knn", "fixed.radius"),
  support.size = NULL,
  radius = NULL,
  min.support = NULL,
  support.buffer = 3L,
  kernel = c("epanechnikov", "triangular", "gaussian", "tricube"),
  coordinate.method = c("coordinates", "local.pca"),
  chart.dim = NULL,
  support.metric = c("auto", "coordinates", "graph.geodesic"),
  auto.chart.support.metric = c("coordinates", "operator", "both"),
  auto.chart.selection.metric = c("coordinates", "operator"),
  row.normalize = c("l2", "none", "l1"),
  local.solver = c("auto", "normal.equations", "qr", "svd"),
  normal.equations.max.condition = 1e+08,
  duplicate.action = c("keep", "error"),
  drop.rank.deficient = TRUE,
  sync.row.normalize = c("l2"),
  sync.min.norm = 1e-12,
  verbose = FALSE,
  ...
)

Arguments

Value

A list of class "slpl_tf_operator" containing A_LPL, C_sync, row metadata, diagnostics, settings, and the embedded "lpl_tf_operator".

Examples

X <- matrix(seq(0, 1, length.out = 18), ncol = 1)
slpl.tf.operator(X, degree = 1L, support.type = "knn",
                 support.size = 7L, kernel = "gaussian")

Extract A Conditional Smoother Matrix

Description

Extract the linear map from responses to predictions for a supported fitted model. R selects the method from the class of object.