testpar: Parametric Monotone/Quadratic vs Smooth Constrained Test

In cgam: Constrained Generalized Additive Model

Parametric Monotone/Quadratic vs Smooth Constrained Test

Description

Performs a hypothesis test comparing a parametric model (e.g., linear, quadratic, warped plane) to a constrained smooth alternative under shape restrictions.

Usage

testpar(
  formula0,
  formula,
  data,
  family = gaussian(link = "identity"),
  ps = NULL,
  edfu = NULL,
  nsim = 200,
  multicore = getOption("cgam.parallel"),
  method = "testpar.fit",
  arp = FALSE,
  p = NULL,
  space = "Q",
  ...
)

Arguments

formula0	A model formula representing the null hypothesis (e.g., linear or quadratic).
formula	A model formula under the alternative hypothesis, possibly with shape constraints.
data	A data frame or environment containing the variables in the model.
family	A description of the error distribution and link function to be used in the model (e.g., gaussian()). Gaussian and binomial are now included.
ps	User-defined penalty term. If NULL, optimal values are estimated.
edfu	User-defined unconstrained effective degrees of freedom.
nsim	Number of simulations to perform to get the mixture distribution of the test statistic. (default is 200).
multicore	Logical. Whether to enable parallel computation (default uses global option cgam.parallel).
method	A character string (default is "testpar.fit"), currently unused.
arp	Logical. If TRUE, uses autoregressive structure estimation.
p	Integer vector or value specifying AR(p) order. Ignored if arp = FALSE.
space	Character. Either "Q" for quantile-based knots or "E" for evenly spaced knots.
...	Additional arguments passed to internal functions.

Value

A list containing test statistics, fitted values, coefficients, and other relevant objects.

Examples

# Example 1: Linear vs. monotone alternative
set.seed(123)
n <- 100
x <- sort(runif(n))
y <- 2 * x + rnorm(n)
dat <- data.frame(y = y, x = x)
ans <- testpar(formula0 = y ~ x, formula = y ~ s.incr(x), multicore = FALSE, nsim = 10)
ans$pval
summary(ans)
## Not run: 
# Example 2: Binary response: linear vs. monotone alternative
set.seed(123)
n <- 100
x <- runif(n)
eta <- 8 * x - 4
mu <- exp(eta) / (1 + exp(eta))
y <- rbinom(n, size = 1, prob = mu)

ans <- testpar(formula0 = y ~ x, formula = y ~ s.incr(x), 
family = binomial(link = "logit"), nsim = 10)
summary(ans)

xs = sort(x)
ord = order(x)
 par(mfrow = c(1,1))
 plot(x, y, cex = .2)
 lines(xs, binomial(link="logit")$linkinv(ans$etahat0)[ord], col=4, lty=4)
 lines(xs, binomial(link="logit")$linkinv(ans$etahat)[ord], col=2, lty=2)
 lines(xs, mu[ord])
 legend("topleft", legend = c("H0 fit (etahat0)", "H1 fit (etahat)", "True mean"),
 col = c(4, 2, 1), lty = c(4, 2, 1), bty = "n", cex = 0.8)

## End(Not run)

# Example 3: Bivariate warped-plane surface test
set.seed(1)
n <- 100
x1 <- runif(n)
x2 <- runif(n)
mu <- 4 * (x1 + x2 - x1 * x2)
eps <- rnorm(n)
y <- mu + eps

# options(cgam.parallel = TRUE) #allow parallel computation
ans <- testpar(formula0 = y ~ x1 * x2, formula = y ~ s.incr.incr(x1, x2), nsim = 10)
summary(ans)

par(mfrow = c(1, 2)) 
persp(ans$k1, ans$k2, ans$etahat0_surf, th=-50, main = 'H0 fit')
persp(ans$k1, ans$k2, ans$etahat_surf, th=-50, main = 'H1 fit')

## Not run: 
# Example 4: AR(1) error, quadratic vs smooth convex
set.seed(123)
n = 100
x = runif(n)
mu = 1+x+2*x^2
eps = arima.sim(n = n, list(ar = c(.4)), sd = 1)
y = mu + eps
ans = testpar(y~x^2, y~s.conv(x), arp=TRUE, p=1, nsim=10)
ans$pval
ans$phi #autoregressive coefficient estimate

## End(Not run)
## Not run: 
# Example 5: Three additive components with shape constraints
n <- 100
x1 <- runif(n)
x2 <- runif(n)
x3 <- runif(n)
z <- rnorm(n)
mu1 <- 3 * x1 + 1
mu2 <- x2 + 3 * x2^2
mu3 <- -x3^2
mu <- mu1 + mu2 + mu3
y <- mu + rnorm(n)

ans <- testpar(formula0 = y ~ x1 + x2^2 + x3^2 + z, 
formula = y ~ s.incr(x1) + s.incr.conv(x2) + s.decr.conc(x3)
 + z, family = "gaussian", nsim = 10)
ans$pval
summary(ans)

## End(Not run)

cgam documentation built on April 12, 2025, 1:09 a.m.