Description Usage Arguments Details See Also Examples
This function tests the equality of nonlinear curves and surface estimations with correlated data based on L2 distance. The semiparametric estimation uses 'gamm4' package with a compond symmetry correlation structure to adjust correlated observations. The specific model considered here is
1 2 | gamm4.grptest(formula, random, test, data, N.boot = 200, m = 225,
parallel = TRUE)
|
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. Note that ids for smooths and fixed smoothing parameters are not supported. |
random |
An optional formula specifying the random effects structure in lmer style. |
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 |
y_ij= m_i(x_ij) + b_i + e_ij,
where m_i(.), are semiparametric smooth functions; b_i are subject-specific random intercept; 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.
gamm4
gam.grptest
plot.gamtest
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 | #Test the equality of three nonlinear curves
m1 <- 120 #number of subjects in group 1
m2 <- 100 #number of subjects in group 2
m3 <- 110 #number of subjects in group 3
n1 <- 3 #number of repeated measurements for each subject in group 1
n2 <- 4 #number of repeated measurements for each subject in group 2
n3 <- 2 #number of repeated measurements for each subject in group 3
sigma1 <- 0.3
sigma2 <- 0.2
sigma3 <- 0.2
sigma.noise1 <- sigma.noise2 <- sigma.noise3 <- 0.1
f1 <- function(u) sin(2*pi*u)
f2 <- f3 <- function(u) sin(2*pi*u)+u/3
N1 <- m1*n1
N2 <- m2*n2
N3 <- m3*n3
x11 <- runif(N1,0,1)
b1i <- rnorm(m1,0,sigma1)
b1 <- rep(b1i, each=n1)
id1 <- rep(1:m1, each=n1)
y1 <- f1(x11) + b1 + rnorm(N1, 0, sigma.noise1)
x21 <- runif(N2,0,1)
b2i <- rnorm(m2,0,sigma2)
b2 <- rep(b2i,each=n2)
id2 <- rep((m1+1):(m1+m2),each=n2)
y2 <- f2(x21) + b2 + rnorm(N2,0,sigma.noise2)
x31 <- runif(N3,0,1)
b3i <- rnorm(m3,0,sigma3)
b3 <- rep(b3i,each=n3)
id3 <- rep((m1+m2+1):(m1+m2+m3),each=n3)
y3 <- f3(x31) + b3 + rnorm(N3,0,sigma.noise2)
dat <- data.frame(rbind(cbind(id1, x11,y1,1), cbind(id2, x21, y2,2), cbind(id3, x31, y3,3)))
colnames(dat)=c('id','x', 'y','grp')
testout <- gamm4.grptest(formula=y~s(x,k=6,bs="cr"), test=~grp,
random=~(1|id), data=dat, N.boot=200, m=225, parallel = TRUE)
testout
plot(testout)
dat0 <- data.frame(rbind(cbind(id3, x31, y3, 3), cbind(id2, x21, y2, 2)))
colnames(dat0)=c('id', 'x', 'y', 'grp')
testout0 <- gamm4.grptest(formula=y~s(x,k=6,bs="cr"), test=~grp,
random=~(1|id), data=dat0, N.boot=200, m=225, parallel= TRUE)
testout0$p.value
plot(testout0, test.statistic = TRUE)
########
## Semiparametric test the equality for regression surfaces with longitudinal data
## Simulate data sets
f1 <- function(u,v) 2*u^2+3*v^2
f2 <- function(u,v) 2*u^2+3*v^2+sin(2*pi*u)
m1 <- 100 #number of subjects in group 1
n1 <- 4 #number of repeated measurements for each subject in group 1
m2 <- 120 #number of subjects in group 2
n2 <- 3 #number of repeated measurements for each subject in group 2
N1 <- m1*n1
N2 <- m2*n2
sigma1 <- 0.2
sigma2 <- 0.15
sigma.noise1 <- 0.04
sigma.noise2 <- 0.05
x11 <- runif(N1,0,1)
x12 <- runif(N1,0,1)
b1i <- rnorm(m1,0,sigma1)
b1 <- rep(b1i,each=n1)
id1 <- rep(1:m1,each=n1)
y1 <- f1(x11,x12) + b1 + rnorm(N1,0, sigma.noise1)
x21 <- runif(N2,0,1)
x22 <- runif(N2,0,1)
b2i <- rnorm(m2,0,sigma2)
b2 <- rep(b2i,each=n2)
id2 <- rep((m1+1):(m1+m2),each=n2)
y2 <- f2(x21,x22) + b2 + rnorm(N2,0,sigma.noise2)
y3 <- f1(x21,x22) + b2 + rnorm(N2,0,sigma.noise2)
dat <- data.frame(rbind(cbind(id1, x11, x12,y1,1), cbind(id2, x21, x22, y2,2)))
colnames(dat)=c('id','x1','x2', 'y','grp')
test.spline1 <- gamm4.grptest(formula=y~t2(x1,x2), test=~grp,
random=~(1|id), data=dat, N.boot=200, m=225, parallel=TRUE)
plot(test.spline1)
plot(test.spline1, type="plotly.persp")
plot(test.spline1, type="plotly.persp", data.pts=TRUE)
dat0 <- data.frame(rbind(cbind(id1, x11, x12 , y1, 1), cbind(id2, x21, x22, y3, 2)))
colnames(dat0)=c('id','x1','x2', 'y','grp')
test.spline0 <- gamm4.grptest(y~t2(x1,x2), test=~grp,
random=~(1|id), data=dat0, N.boot=200, m=225, parallel=TRUE)
test.spline0
plot(test.spline0, test.statistic = FALSE)
plot(test.spline0)
plot(test.spline0, type="plotly.persp")
########
## Data analyses with internal "outchild" dataset
data("outchild")
outchild1016 <- outchild[(outchild$age<=16 & outchild$age>10),]
child.repw <- outchild1016[(outchild1016$RACE==1),]
child.reptest1 <- gamm4.grptest(HEIGHT~s(age), random=~(1|SID),
test=~SEX, data=child.repw, parallel = TRUE)
child.reptest1
plot(child.reptest1)
plot(child.reptest1,test.statistic = FALSE)
child.reptest2 <- gamm4.grptest(WEIGHT~t2(age,HEIGHT), random=~(1|SID),
test=~SEX, data = child.repw, parallel = TRUE)
plot(child.reptest2,type="plotly.persp")
plot(child.reptest2,type="contour")
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.