pif.confidence: Confidence Intervals for the Potential Impact Fraction

In pifpaf: Potential Impact Fraction and Population Attributable Fraction for Cross-Sectional Data

Description

Function that estimates confidence intervals for the Potential Impact Fraction pif from a cross-sectional sample of the exposure X with a known Relative Risk function rr with parameter theta, where the Potential Impact Fraction is given by:

PIF = (mean(rr(X; theta)) - mean(rr(cft(X);theta)))/mean(rr(X; theta)) .

Usage

pif.confidence(X, thetahat, rr, thetavar, cft = NA, method = "empirical",
  confidence_method = "bootstrap", confidence = 95, nsim = 1000,
  weights = rep(1/nrow(as.matrix(X)), nrow(as.matrix(X))), Xvar = var(X),
  deriv.method.args = list(), deriv.method = "Richardson", adjust = 1,
  n = 512, ktype = "gaussian", bw = "SJ", check_exposure = TRUE,
  check_cft = TRUE, check_rr = TRUE, check_xvar = TRUE,
  check_integrals = TRUE, check_thetas = TRUE, is_paf = FALSE)

Arguments

X	Random sample (data.frame) which includes exposure and covariates or sample mean if "approximate" method is selected.
thetahat	Asymptotically consistent or Fisher consistent estimator (vector) of theta for the Relative Risk function. thetahat should be asymptotically normal (N(theta, var_of_theta)) with mean theta and estimated variance var_of_theta.
rr	function for Relative Risk which uses parameter theta. The order of the parameters shound be rr(X, theta).
thetavar	Estimator of variance var_of_theta. **Optional**
cft	Function cft(X) for counterfactual. Leave empty for the Population Attributable Fraction paf where counterfactual is that of the theoretical minimum risk exposure X0 such that rr(X0,theta) = 1.
method	Either "empirical" (default), "kernel" or "approximate". For details on estimation methods see pif.
confidence_method	Either bootstrap (default), linear, loglinear. See paf details for additional information.
confidence	Confidence level % (default 95).
nsim	Number of simulations for estimation of variance.
weights	Normalized survey weights for the sample X.
Xvar	Variance of exposure levels (for "approximate" method).
deriv.method.args	method.args for hessian (for "approximate" method).
deriv.method	method for hessian. Don't change this unless you know what you are doing (for "approximate" method).
adjust	Adjust bandwith parameter (for "kernel" method) from density.
n	Number of equally spaced points at which the density (for "kernel" method) is to be estimated (see density).
ktype	kernel type: "gaussian", "epanechnikov", "rectangular", "triangular", "biweight", "cosine", "optcosine" (for "kernel" method). Additional information on kernels in density.
bw	Smoothing bandwith parameter (for "kernel" method) from density. Default "SJ".
check_exposure	boolean Check that exposure X is positive and numeric.
check_cft	boolean Check that counterfactual function cft reduces exposure.
check_rr	boolean Check that Relative Risk function rr equals 1 when evaluated at 0.
check_xvar	boolean Check Xvar is covariance matrix.
check_integrals	boolean Check that counterfactual cft and relative risk's rr expected values are well defined for this scenario.
check_thetas	boolean Check that theta associated parameters are correctly inputed for the model.
is_paf	Boolean forcing evaluation of paf. This forces the pif function to ignore the inputed counterfactual and set it to the theoretical minimum risk value of rr = 1.

Value

pifvec Vector with lower ("Lower_CI"), and upper ("Upper_CI") confidence bounds for the pif as well as point estimate "Point_Estimate" and estimated variance of log(pif) (if confidence_method is "loglinear").

Note

For more information on kernels see density.

Do not use the $ operator when using "approximate" method.

Author(s)

Rodrigo Zepeda-Tello [e-mail]

Dalia Camacho-Garc<c3><ad>a-Forment<c3><ad> [e-mail]

See Also

See paf.confidence for confidence interval estimation of paf, and pif for only point estimate.

Sensitivity analysis plots can be done with paf.plot, and paf.sensitivity

Examples

#Example 1: Exponential Relative Risk
#--------------------------------------------
set.seed(18427)
X        <- data.frame(Exposure = rnorm(100,3,1))
thetahat <- 0.12
thetavar <- 0.02
rr       <- function(X, theta){exp(theta*X)}


#Counterfactual of halving exposure
cft   <- function(X){ 0.5*X }

#Using bootstrap method
pif.confidence(X, thetahat, rr, thetavar, cft)

## Not run: 
#Same example with loglinear method
pif.confidence(X, thetahat, rr, thetavar, cft, confidence_method = "loglinear")

#Same example with linear method (usually the widest and least precise)
pif.confidence(X, thetahat, rr, thetavar, cft, confidence_method = "linear")


#Example 2: Linear Relative Risk
#--------------------------------------------
set.seed(18427)
X        <- data.frame(Exposure = rbeta(100,3,1))
thetahat <- 0.17
thetavar <- 0.01
rr       <- function(X, theta){theta*X + 1}
cft      <- function(X){ 0.5*X }
weights             <- runif(100)
normalized_weights  <- weights/sum(weights)
pif.confidence(X, thetahat, rr, thetavar, cft, weights = normalized_weights)

#Same example with more complex counterfactual that reduces 
#only the values > 0.75 are halved
cft       <- function(X){

   #Indentify the ones with "a lot" of exposure:
   where_excess_exposure    <- which(X[,"Exposure"] > 0.75)             
   
   #Halve their exposure
   X[where_excess_exposure, "Exposure"] <- 
            X[where_excess_exposure, "Exposure"]/2  
   return(X)
}
pif.confidence(X, thetahat, rr, thetavar, cft, weights = normalized_weights)


#Example 3: Multivariate Linear Relative Risk
#--------------------------------------------
set.seed(18427)
X1       <- rnorm(100,4,1)
X2       <- rnorm(100,2,0.4)
thetahat <- c(0.12, 0.03)
thetavar <- diag(c(0.01, 0.02))

#But the approximate method crashes due to operator
Xmean <- data.frame(Exposure = mean(X1), 
                    Covariate = mean(X2))
Xvar  <- var(cbind(X1, X2))

#When creating relative risks avoid using the $ operator
#as it doesn't work under approximate method of PAF
rr_not    <- function(X, theta){
               exp(theta[1]*X$Exposure + theta[2]*X$Covariate)
             }
rr_better <- function(X, theta){
               exp(theta[1]*X[,"Exposure"] + theta[2]*X[,"Covariate"])
             }
             
#Creating a counterfactual. 
cft  <- function(X){
   Y               <- X
   Y[,"Exposure"]  <- 0.5*X[,"Exposure"]
   Y[,"Covariate"] <- 1.1*X[,"Covariate"] + 1
   return(Y)
}

pif.confidence(Xmean, thetahat, rr_better, thetavar, cft, 
method = "approximate", Xvar = Xvar) 

## End(Not run)

## Not run: 
#Error: $ operator in rr definitions don't work in approximate
pif.confidence(Xmean, thetahat, rr_not, thetavar, cft, method = "approximate", Xvar = Xvar)

## End(Not run)

## Not run: 
#Example 4: Categorical Relative Risk & Exposure
#--------------------------------------------
set.seed(18427)
mysample  <- sample(c("Normal","Overweight","Obese"), 100, 
                   replace = TRUE, prob = c(0.4, 0.1, 0.5))
X        <- data.frame(Exposure = mysample)

thetahat <- c(1, 1.2, 1.5)
thetavar <- diag(c(0.1, 0.2, 0.3))

#Categorical relative risk function
rr <- function(X, theta){

   #Create return vector with default risk of 1
   r_risk <- rep(1, nrow(X))
   
   #Assign categorical relative risk
   r_risk[which(X[,"Exposure"] == "Normal")]      <- thetahat[1]
   r_risk[which(X[,"Exposure"] == "Overweight")]  <- thetahat[2]
   r_risk[which(X[,"Exposure"] == "Obese")]       <- thetahat[3]
   
   return(r_risk)
}


#Counterfactual of reducing all obesity to normality
cft <- function(X){
   X[which(X[,"Exposure"] == "Obese"),] <- "Normal"
   return(X)
}

pif.confidence(X, thetahat, rr, thetavar, cft, check_rr = FALSE)


#Example 5: Categorical Relative Risk & continuous exposure
#----------------------------------------------------------
set.seed(18427)
BMI      <- data.frame(Exposure = rlnorm(100, 3.1, sdlog = 0.1))

#Theoretical minimum risk exposure is 20kg/m^2 in borderline "Normal" category
BMI_adjusted <- BMI - 20

thetahat <- c(Malnourished = 2.2, Normal = 1, Overweight = 1.8, 
              Obese = 2.5)
thetavar <- diag(c(0.1, 0.2, 0.2, 0.1))

rr       <- function(X, theta){
     
     #Create return vector with default risk of 1
     r_risk <- rep(1, nrow(X))
   
     #Assign categorical relative risk
     r_risk[which(X[,"Exposure"] < 0)]             <- theta[1] #Malnourished
     r_risk[intersect(which(X[,"Exposure"] >= 0), 
                      which(X[,"Exposure"] < 5))]  <- theta[2] #Normal
     r_risk[intersect(which(X[,"Exposure"] >= 5), 
                      which(X[,"Exposure"] < 10))] <- theta[3] #Overweight
     r_risk[which(X[,"Exposure"] >= 10)]           <- theta[4] #Obese
   
   return(r_risk)
}

#Counterfactual of everyone in normal range
cft <- function(bmi){
     bmi           <- data.frame(rep(2.5, nrow(bmi)), ncol = 1)
     colnames(bmi) <- c("Exposure")
     return(bmi)
}

pif.confidence(BMI_adjusted, thetahat, rr, thetavar, cft, 
                check_exposure = FALSE, method = "empirical")


#Example 6: Bivariate exposure and rr ("classical PAF")
#------------------------------------------------------------------
set.seed(18427)
mysample  <- sample(c("Exposed","Unexposed"), 1000, 
                replace = TRUE, prob = c(0.1, 0.9))
X         <- data.frame(Exposure = mysample)
theta     <- c("Exposed" = 2.5, "Unexposed" = 1.2)   
thetavar  <- matrix(c(0.04, 0.02, 0.02, 0.03), ncol = 2)
rr        <- function(X, theta){
   
   #Create relative risk function
   r_risk <- rep(1, nrow(X))
   
   #Assign values of relative risk
   r_risk[which(X[,"Exposure"] == "Unexposed")] <- theta["Unexposed"]
   r_risk[which(X[,"Exposure"] == "Exposed")]   <- theta["Exposed"]
   
   return(r_risk)
}      

#Counterfactual of reducing the exposure in half of the individuals
cft <- function(X){

   #Find out which ones are exposed
   Xexp  <- which(X[,"Exposure"] == "Exposed")
   
   #Use only half of the exposed randomly
   reduc <- sample(Xexp, length(Xexp)/2)
   
   #Unexpose those individuals
   X[reduc, "Exposure"] <- "Unexposed"
   
   return(X)
}

pif.confidence(X, theta, rr, thetavar, cft)

#Example 7: Continuous exposure, several covariates
#------------------------------------------------------------------
X <- data.frame(Exposure = rbeta(100, 2, 3),
                Age      = runif(100, 20, 100),
                Sex      = sample(c("M","F"), 100, replace = TRUE),
                BMI      = rlnorm(100, 3.2, 0.2))
thetahat <- c(-0.1, 0.05, 0.2, -0.4, 0.3, 0.1)

#Create variance of theta
almostvar <- matrix(runif(6^2), ncol = 6)

thetavar <- t(almostvar) %*% almostvar

rr <- function(X, theta){
     #Create risk vector
     Risk    <- rep(1, nrow(X))
     
     #Identify subpopulations
     males   <- which(X[,"Sex"] == "M")
     females <- which(X[,"Sex"] == "F")
     
     #Calculate population specific rr
     Risk[males] <- theta[1]*X[males,"Exposure"] + 
                                      theta[2]*X[males,"Age"]^2 + 
                                      theta[3]*X[males,"BMI"]/2 
                                     
     Risk[females] <- theta[4]*X[females,"Exposure"] + 
                                      theta[5]*X[females,"Age"]^2 + 
                                      theta[6]*X[females,"BMI"]/2 
                                     
    return(Risk)
}

#Counterfactual of reducing BMI
cft <- function(X){
    excess_bmi           <- which(X[,"BMI"] > 25)
    X[excess_bmi,"BMI"]  <- 25
    return(X)
}

pif.confidence(X, thetahat, rr, thetavar, cft)

## End(Not run)

pifpaf documentation built on May 1, 2019, 9:11 p.m.