dot-lagrangian_project: Lagrangian Projection via OLS Reformulation

In lgspline: Lagrangian Multiplier Smoothing Splines for Smooth Function Estimation

Lagrangian Projection via OLS Reformulation

Description

Given unconstrained (or GEE-initialized) coefficient estimates \hat{\boldsymbol{\beta}}, computes the constrained estimate \tilde{\boldsymbol{\beta}} = \mathbf{U}\hat{\boldsymbol{\beta}} where \mathbf{U} = \mathbf{I} - \mathbf{G}\mathbf{A} (\mathbf{A}^{\top}\mathbf{G}\mathbf{A})^{-1}\mathbf{A}^{\top}.

Usage

.lagrangian_project(
  GhalfXy,
  Ghalf,
  A,
  K,
  p_expansions,
  R_constraints,
  constraint_value_vectors,
  family,
  parallel_aga,
  parallel_matmult,
  cl,
  chunk_size,
  num_chunks,
  rem_chunks,
  parallel_qr = FALSE
)

Arguments

GhalfXy	Numeric column vector \mathbf{y}^* of length P; produced by multiplying \mathbf{G}^{1/2} (Path 2) or \mathbf{G}^{-1/2} (Path 3) into the unconstrained estimate.
Ghalf	List of \mathbf{G}^{1/2}_k matrices (one per partition).
A	Constraint matrix \mathbf{A} (P \times R).
K	Integer; number of interior knots.
p_expansions	Integer; number of basis terms per partition.
R_constraints	Integer; number of columns of \mathbf{A}.
constraint_value_vectors	List of constraint right-hand-side vectors encoding \mathbf{A}^{\top}\boldsymbol{\beta} = \mathbf{c}. When non-empty and nonzero, the particular solution is added.
family	GLM family object (used for linkinv(0) fallback when NA values arise in \mathbf{y}^*).
parallel_aga, parallel_matmult, parallel_qr	Logical flags for parallel computation.
cl	Parallel cluster object.
chunk_size, num_chunks, rem_chunks	Parallel distribution parameters.

Details

Rather than forming \mathbf{U} directly, the projection is reformulated as a residual from an OLS problem (the \mathbf{G}^{1/2}\mathbf{r}^* trick):

\mathbf{y}^* = \mathbf{G}^{-1/2}\hat{\boldsymbol{\beta}}

\mathbf{X}^* = \mathbf{G}^{1/2}\mathbf{A}

\mathbf{r}^* = (\mathbf{I} - \mathbf{X}^*(\mathbf{X}^{*\top} \mathbf{X}^*)^{-1}\mathbf{X}^{*\top})\mathbf{y}^*

\tilde{\boldsymbol{\beta}} = \mathbf{G}^{1/2}\mathbf{r}^*

This avoids explicitly forming and inverting \mathbf{A}^{\top}\mathbf{G}\mathbf{A}; the most expensive step is the QR decomposition of the R \times R system inside .lm.fit, which is far cheaper than the full P \times P solve.

Value

A list of length K+1, each element a column vector of constrained coefficients \tilde{\boldsymbol{\beta}}_k.

lgspline documentation built on May 8, 2026, 5:07 p.m.