gamm4.test (version 0.1.0)

gam.grptest: Test the equality of nonlinear curves and surface estimations by semiparametric method

Description

This function tests the equality of nonlinear curves and surface estimations based on L2 distance. The semiparametric estimation uses 'mgcv' package. The specific model considered here is

Usage

gam.grptest(formula, test, data, N.boot = 200, m = 225, parallel = FALSE)

Arguments

formula

A GAM formula. This is like the formula for a glm except that smooth terms (s and t2 but not te) can be added to the right hand side of the formula.

test

An indicator of variable for testing nonlinear curves or surface estimations

data

A data frame or list containing the model response variable and covariates required by the formula.

N.boot

the number of bootstrap replicates. This should be a single positive integer.

m

the number of the sampling points for the Monte-Carlo integration.

parallel

Parallel computation of semiparametric estimations with bootstrap samples for getting test statistics under null hypothesis.

Details

y_ij= m_i(x_ij) + e_ij,

where m_i(.), are semiparametric smooth functions; e_ij are subject-specific errors. The errors e_ij do not have to be independent N(0, sigma^2) errors. The errors can be heteroscedastic, i.e., e_ij = sigma_i(x_ij) * u_ij, where u_ij are independent identically distributed errors with mean 0 and variance 1.

We are interested in the problem of testing the equality of the regression curves (when x is one-dimensional) or surfaces (when x is two-dimensional),

H_0: m_1(.) = m_2(.) = ... v.s. H_1: otherwise

The problem can also be viewed as the test of the equality in the one-sample problem for functional data.

A bootstrap algorithm is applied to test the equality of semiparametric curves or surfaces based on L2 distance.

See Also

gam gamm4.grptest plot.gamtest T.L2c

Examples

# NOT RUN {
# }
# NOT RUN {
n1 <- 200
x1 <- runif(n1,min=0, max=3)
sd1 <- 0.2
e1 <- rnorm(n1,sd=sd1)
y1 <- sin(2*x1) + cos(2*x1) + e1

n2 <- 120
x2 <- runif(n2, min=0, max=3)
sd2 <- 0.25
e2 <- rnorm(n2, sd=sd2)
y2 <- sin(2*x2) + cos(2*x2) + x2 + e2

data.bind <- rbind(cbind(x1,y1,1), cbind(x2,y2,2))
data.bind <- data.frame(data.bind)
colnames(data.bind)=c('x','y','group')

t1 <- gam.grptest(y~s(x,bs="cr"), test=~group, data=data.bind, parallel=FALSE)
t1
plot(t1)

########
## Semiparametric test the equality for regression surfaces
## Simulate data sets

n1 <- 500
x11 <- runif(n1,min=0, max=3)
x12 <- runif(n1,min=0, max=3)
sd1 <- 0.2
e1 <- rnorm(n1,sd=sd1)
y1 <- 2*x11^2 + 3*x12^2  + e1

n2 <- 420
x21 <- runif(n2, min=0, max=3)
x22 <- runif(n2, min=0, max=3)
sd2 <- 0.25
e2 <- rnorm(n2, sd=sd2)
y2 <- 2*x21^2 + 3*x22^2 + 6*sin(2*pi*x21) + e2

n3 <- 550
x31 <- runif(n3,min=0, max=3)
x32 <- runif(n3,min=0, max=3)
sd3 <- 0.2
e3 <- rnorm(n3,sd=sd1)
y3 <- 2*x31^2 + 3*x32^2  + e3

data.bind <- rbind(cbind(x11, x12 ,y1,1), cbind(x21, x22, y2,2), cbind(x31, x32, y3,3))
data.bind <- data.frame(data.bind)
colnames(data.bind)=c('x1','x2', 'y','group')

tspl <- gam.grptest(y~s(x1,x2), test=~group, data=data.bind, N.boot=200, m=225, parallel=FALSE)
tspl$p.value #p-value
plot(tspl, test.statistic = TRUE)
plot(tspl, type="contour")
plot(tspl, type="persp")
plot(tspl, type="plotly.persp")
plot(tspl, type="plotly.persp",data.pts=TRUE)

########
## Data analyses with internal "outchild" dataset

data("outchild")
child<- outchild[order(outchild$SID,outchild$age),]
bs <- aggregate(.~SID, child, FUN=head, 1)

childcur <- bs[,c("SEX","WEIGHT","age")]
test.grpsex1 <- gam.grptest(WEIGHT~s(age), test=~SEX, data=childcur)
test.grpsex1
plot(test.grpsex1)
plot(test.grpsex1,test.statistic=TRUE)

childsurf <- bs[,c("SEX","HEIGHT","WEIGHT","age")]
test.grpsex2 <- gam.grptest(WEIGHT~s(HEIGHT,age), test=~SEX, data=childsurf)
test.grpsex2
plot(test.grpsex2)
plot(test.grpsex2, type="plotly.persp")
plot(test.grpsex2, type="plotly.persp",data.pts=TRUE)
# }