Description Usage Arguments Value Examples
This is the base function for testing average treatment effects. Users can use specify the nuisance function values by themselves.
1 2 3 4 5 6 7 8 9 10 11 12 13 | drdrtest.base(
y,
a,
pi,
varpi,
mu,
ma,
arange,
h = NULL,
b = 1000,
dist = "TwoPoint",
a.grid.size = 401
)
|
y |
A vector containing the outcomes for each observation |
a |
A vector containing the treatment levels (dosage) for each observation |
pi |
A vector containing the propensity scores for each observation |
varpi |
A vector containing the mean propensity scores for each observation |
mu |
A vector containing the outcome regression function values for each observation |
ma |
A vector containing the mean outcome regression fucntion values for each observation |
arange |
A vector of length 2 giving the lower bound and upper bound of treatment levels |
h |
bandwidth to be used in kernel regression. If not specified, will by default use "rule of thumb" bandwidth selector |
b |
number of Bootstrap samples to be generated |
dist |
distibution used to generate residuals for Bootstrap samples. Currently only have two options, "TwoPoint" and "Rademachar" |
a.grid.size |
size of equally spaced grid points over |
A list containing
P value of the test result
Value of the observed test statistic
A vector containing test statistic values from Bootstrap samples
A list containg evalution points of average treatment effect and the corresponding values
Bandwidth used in kernel regression
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 | mu.mod<-function(a,l,delta,height){
mu <- as.numeric(l%*%c(0.2,0.2,0.3,-0.1))+triangle(a-2.5,delta,height)+a*(-0.1*l[,1]+0.1*l[,3])
return(mu)
}
triangle <- function(a,delta,height){
y <- exp(-a^2/((delta/2)^2))*height
return(y)
}
set.seed(2000)
n <- 500
d <- 4
sigma <- 0.5
delta <- 1
height <- 0
arange<-c(0.01,4.99)
l <- matrix(rnorm(n*d),ncol=d)
colnames(l) <- paste("l",1:4,sep="")
logit.lambda <- as.numeric(l%*%c(0.1,0.1,-0.1,0.2))
lambda <- exp(logit.lambda)/(1+exp(logit.lambda))
a <- rbeta(n, shape1 = lambda, shape2 =1-lambda)*5
mu <- mu.mod(a,l,delta,height)
residual.list <- rnorm(n,mean=0,sd=sigma)
y <- mu+residual.list
## We use the oracal propensity score and outcome regression for illustration
pilist <- dbeta(a/5, shape1=lambda, shape2 = 1-lambda)/5
varpilist <- colMeans(matrix(dbeta(rep(a,each=n)/5,
shape1=rep(lambda,n),
shape2 = 1-rep(lambda,n))/5, nrow=n))
mulist <- mu
malist <-colMeans(matrix(mu.mod(rep(a,each=n),l[rep(1:n,n),],delta,height),nrow=n))
out <- drdrtest.base(y,a,pilist,varpilist,mulist,malist,arange)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.