Description Usage Arguments Value Examples
This is the function for testing average treatment effects with user specified nuisance functions.
1 2 3 4 5 6 7 8 9 10 11 12 13 |
y |
A vector containing the outcomes for each observation |
a |
A vector containing the treatment levels (dosage) for each observation |
l |
A data.frame containing the observations of covariates |
arange |
A vector of length 2 giving the lower bound and upper bound of treatment levels |
pifunc |
A user specifid function or wapper that takes treatment a as the first argument and covariates l as the second argument and return propensit scores |
mufunc |
A user specifid function or wapper that takes treatment a as the first argument and covariates l as the second argument and return outcome regression values |
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" |
pi.low |
Lower bound to truncate propensity scores |
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 36 37 38 39 | 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.05
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
pifunc <- function(a,l){
l <- as.matrix(l)
logit.lambda <- as.numeric(l%*%c(0.1,0.1,-0.1,0.2))
lambda <- exp(logit.lambda)/(1+exp(logit.lambda))
return(dbeta(a/5,shape1=lambda,shape2 = 1-lambda)/5)
}
mufunc <- function(a,l){
l <- as.matrix(l)
return(mu.mod(a,l,delta,height))
}
out <- drdrtest(y,a,data.frame(l),arange,pifunc,mufunc)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.