ivtools (version 2.3.0)

ivcoxph: Instrumental variable estimation of the causal exposure effect in Cox proportional hazards (PH) models

Description

ivcoxph performs instrumental variable estimation of the causal exposure effect in Cox PH models with individual-level data. Below, \(Z\), \(X\), and \(T\) are the instrument, the exposure, and the outcome, respectively. \(L\) is a vector of covariates that we wish to control for in the analysis; these would typically be confounders for the instrument and the outcome.

Usage

ivcoxph(estmethod, X, fitZ.L=NULL, fitX.LZ=NULL, fitX.L=NULL, fitT.LX=NULL, 
  fitT.LZX=NULL, data, formula=~1, ctrl=FALSE, clusterid=NULL, t=NULL, 
  vcov.fit=TRUE, ...)

Arguments

estmethod

a string specifying the desired estimation method; either "ts" for two-stage estimation, or "g" for G-estimation.

X

a string specifying the name of the exposure \(X\) in data. This is not needed if fitX.LZ is specified.

fitZ.L

an object of class "glm", as returned by the glm function in the stats package. This is a fitted GLM for \(E(Z|L)\). If there are no covariates, then fitZ.L may be specified as a model with an intercept only. This argument is not used when estmethod="ts".

fitX.LZ

an object of class "glm", as returned by the glm function in the stats package. This is a fitted GLM for \(E(X|L,Z)\).

fitX.L

an object of class "glm", as returned by the glm function in the stats package. This is a fitted GLM for \(E(X|L)\). If there are no covariates, then fitX.L may be specified as a model with an intercept only. This argument is not used when estmethod="ts".

fitT.LX

an object of class "coxph", as returned by the coxph function in the survival package. This is a fitted Cox PH model for \(\lambda(t|L,X)\). This argument is not used when estmethod="g".

fitT.LZX

either an object of class "coxph" or an object of class "survfit", as returned by the coxph function in the survival package. This is a fitted Cox PH model for \(\lambda(t|L,Z,X)\) or a non-parametric model for \(S(t|L,Z,X)\), respectively. This argument is not used when estmethod="ts".

data

a data frame containing the variables in the model. The covariates, instrument, exposure and outcome can have arbitrary names, e.g. they don't need to be called L, Z, X and T.

formula

an object of class "formula", with no left-hand side. This specifies the causal interaction terms \(m(L)\); see `Details'. Defaults to ~1, i.e. main effect only. This argument is not used when estmethod="ts".

ctrl

logical. Should the control function \(R=X-\hat{X}\) be used when re-fitting fitT.LX? This argument is not used when estmethod="g".

clusterid

an optional string containing the name of a cluster identification variable when data are clustered. Specifying clusterid corrects the standard errors but does not affect the estimates.

t

a numeric scalar specifying the time point at which to solve the estimating equation when estmethod="g"; see `Details'. If not specified, then the estimating equation is solved at the optimal value of t, defined as the value that minimizes \(trace\{var(\hat{\psi})\}\); see Martinussen et al (2017). This argument is not used when estmethod="ts".

vcov.fit

logical. Should the variance-covariance matrix be computed?

...

optional arguments passed on to the nleqslv function, which is used to solve the estimating equations when estmethod="g". See the help pages for nleqslv. This argument is not used when estmethod="ts".

Value

ivcoxph returns an object of class "ivcoxph", which inherits from class "ivmod". An object of class "ivcoxph" is a list containing

call

the matched call.

input

input is a list containing all input arguments

est

a vector containing the estimate of \(\psi\).

vcov

the variance-covariance matrix for the estimate of \(\psi\), obtained with the sandwich formula.

estfunall

a matrix of all subject-specific contributions to the estimating functions used in the estimation process. One row for each subject, one column for each parameter. If estmethod="ts", then the first columns correspond to the parameters estimated by fitX.LZ, and the last columns correspond to the parameters estimated by the re-fitted model fitT.LX. If estmethod="g", then the first columns correspond to \(\psi\), and the remaining columns correspond the parameters estimated by fitZ.L, fitX.LZ, fitX.L and fitT.LZX, whichever were used in the estimation process.

d.estfun

the jacobian matrix of colMeans(estfunall).

converged

logical. Was a solution found to the estimating equations?

fitT.LX

the re-fitted model fitT.LX used in the estimation process when estmethod="ts". This element is NULL when estmethod="g".

t

the value of t used in the estimation process. This element is NULL when estmethod="ts".

Details

ivcoxph estimates the parameter \(\psi\) in the causal Cox PH model $$\textrm{log}\{\lambda(t|L,Z,X)\}-\textrm{log}\{\lambda_0(t|L,Z,X)\}=m^T(L)X\psi.$$ Here, \(\lambda_0(t|L,Z,X)\) is counterfactual hazard function, had the exposure been set to 0. The vector function \(m(L)\) contains interaction terms between \(L\) and \(X\). If estmethod="ts", then these are specified implicitly through the model fitT.LX. If estmethod="g", then these are specified explicitly through the formula argument.

If estmethod="ts", then two-stage estimation of \(\psi\) is performed. In this case, the model fitX.LZ is used to construct predictions \(\hat{X}=\hat{E}(X|L,Z)\). These predictions are subsequently used to re-fit the model fitT.LX, with \(X\) replaced with \(\hat{X}\). The obtained coefficient(s) for \(\hat{X}\) in the re-fitted model is the two-stage estimator of \(\psi\).

If estmethod="g", then G-estimation of \(\psi\) is performed. In this case, the estimator is obtained as the solution to the estimating equation $$H(\psi)=\sum_{i=1}^n\hat{d}(L_i,Z_i)h_i(\psi;t)=0,$$ where $$h_i(\psi;t)=\hat{S}(t|L_i,Z_i,X_i)^{\textrm{exp}\{-m^T(L_i)\psi X_i\}}.$$ The estimated function \(\hat{d}(L,Z)\) is chosen so that the true function has conditional mean 0, given \(L\); \(E\{d(L,Z)|L\}=0\). The specific form of \(\hat{d}(L,Z)\) is determined by the user-specified models. If fitX.LZ and fitX.L are specified, then \(\hat{d}(L,Z)=m(L)\{\hat{E}(X|L,Z)-\hat{E}(X|L)\}\), where \(\hat{E}(X|L,Z)\) and \(\hat{E}(X|L)\) are obtained from fitX.LZ and fitX.L, respectively. If these are not specified, then \(\hat{d}(L,Z)=m(L)\{Z-\hat{E}(Z|L)\}\), where \(\hat{E}(Z|L)\) is obtained from fitZ.L, which then must be specified. The estimating equation is solved at the value of \(t\) specified by the argument t. \(\hat{S}(t|L_i,Z_i,X_i)\) is an estimate of \(S(t|L_i,Z_i,X_i)\) obtained from the model fitT.LZX.

References

Martinussen T., Sorensen D.D., Vansteelandt S. (2019). Instrumental variables estimation under a structural Cox model. Biostatistics 20(1), 65-79.

Sjolander A., Martinussen T. (2019). Instrumental variable estimation with the R package ivtools. Epidemiologic Methods 8(1), 1-20.

Tchetgen Tchetgen E.J., Walter S., Vansteelandt S., Martinussen T., Glymour M. (2015). Instrumental variable estimation in a survival context. Epidemiology 26(3), 402-410.

Examples

# NOT RUN {
require(survival)

set.seed(9)

##Note: the parameter values in the examples below are chosen to make 
##Y0 independent of Z, which is necessary for Z to be a valid instrument.

n <- 10000
psi0 <- 0.5
Z <- rbinom(n, 1, 0.5)
X <- rbinom(n, 1, 0.7*Z+0.2*(1-Z))
m0 <- exp(0.8*X-0.41*Z) #T0 independent of Z at t=1
T <- rexp(n, rate=exp(psi0*X+log(m0)))
C <- rexp(n, rate=exp(psi0*X+log(m0))) #50% censoring
d <- as.numeric(T<C)
T <- pmin(T, C)
data <- data.frame(Z, X, T, d)

#two-stage estimation
fitX.LZ <- glm(formula=X~Z, data=data)
fitT.LX <- coxph(formula=Surv(T, d)~X, data=data)
fitIV <- ivcoxph(estmethod="ts", fitX.LZ=fitX.LZ, fitT.LX=fitT.LX, data=data, 
  ctrl=TRUE)
summary(fitIV)

#G-estimation with non-parametric model for S(t|L,Z,X) and model for Z
fitZ.L <- glm(formula=Z~1, data=data)
fitT.LZX <- survfit(formula=Surv(T, d)~X+Z, data=data)
fitIV <- ivcoxph(estmethod="g", X="X", fitZ.L=fitZ.L, fitT.LZX=fitT.LZX, 
  data=data, t=1)
summary(fitIV)

#G-estimation with Cox model for \lambda(t|L,Z,X) and model for Z
fitZ.L <- glm(formula=Z~1, data=data)
fitT.LZX <- coxph(formula=Surv(T, d)~X+X+X*Z, data=data)
fitIV <- ivcoxph(estmethod="g", X="X", fitZ.L=fitZ.L, fitT.LZX=fitT.LZX, 
  data=data, t=1)
summary(fitIV)

# }