Learn R Programming

cgam (version 1.24)

testpar: 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",
  ...
)

Value

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

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.

Examples

Run this code

# 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)
if (FALSE) {
# 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)
}

# 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')

if (FALSE) {
# 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
}
if (FALSE) {
# 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)
}

Run the code above in your browser using DataLab