PK
nTBWíöíö
MR_survival_sim_rev2.texUT ñœeñœe%\documentclass[12pt]{article}
\documentclass[AMA,STIX1COL]{WileyNJD-v2}
\articletype{Article Type}%
\received{}
\revised{}
\accepted{}
\raggedbottom
\usepackage{graphicx, amsmath, amsfonts, multirow, textcomp, lscape, subfigure, color, float, centernot, setspace}
%\usepackage[left=2.5cm,top=3.5cm, bottom=3.5cm,right=2.5cm]{geometry}
%\usepackage[colorlinks,citecolor=blue,urlcolor=blue]{hyperref}
%\usepackage{amsmath}
%\usepackage{enumerate}
%\usepackage[round]{natbib}
%\usepackage{url} % not crucial - just used below for the URL
%\renewcommand{\abstractname}{Summary}
\newcommand{\bs}[1]{\boldsymbol{#1}}
\newcommand{\mb}[1]{\mathbf{#1}}
\newcommand{\bx}{\mathbb{X}}
\newcommand{\wt}[1]{\widetilde{#1}}
\newcommand{\rc}{\color{red}}
\newcommand{\bc}{\color{blue}}
\renewcommand{\theenumi}{(A\arabic{enumi})}
\renewcommand{\labelenumi}{\theenumi}
\newcommand\indep{\protect\mathpalette{\protect\independenT}{\perp}}
\def\independenT#1#2{\mathrel{\rlap{$#1#2$}\mkern2mu{#1#2}}}
\newcommand{\notindep}{\centernot{\indep}}
\def\E{\text{E}}
\def\Cov{\text{Cov}}
\def\bA{\mb{A}}
\def\bB{\mb{B}}
\def\bb{\mb{b}}
\def\ba{\mb{a}}
\def\bC{\mb{C}}
\def\bD{\mb{D}}
\def\bG{\mb{G}}
\def\bH{\mb{H}}
\def\bI{\mb{I}}
\def\bK{\mb{K}}
\def\bM{\mb{M}}
\def\bS{\mb{S}}
\def\bU{\mb{U}}
\def\bV{\mb{V}}
\def\bW{\mb{W}}
\def\bX{\mb{X}}
\def\bZ{\mb{Z}}
\def\bXi{\mb{X}_i}
\def\bDi{\mb{D}_i}
\def\bSi{\mb{S}_i}
\def\bZi{\mb{Z}_i}
\def\bWi{\mb{W}_i}
\def\mS{\mathcal{S}}
\def\mG{\mathcal{G}}
\def\q3uad{\hskip3em\relax}
\def\q10uad{\hskip10em\relax}
\newcommand{\balpha}{\bs{\alpha}}
\newcommand{\bbeta}{\bs{\beta}}
\newcommand{\bgamma}{\bs{\gamma}}
\newcommand{\bLambda}{\bs{\Lambda}}
\newcommand{\bDelta}{\bs{\Delta}}
\newcommand{\bOmega}{\bs{\Omega}}
\newcommand{\blambda}{\bs{\lambda}}
\newcommand{\btheta}{\bs{\theta}}
\newcommand{\bzeta}{\bs{\zeta}}
\newcommand{\bpsi}{\bs{\psi}}
%\theoremstyle{plain}
%\newtheorem{theorem}{Theorem}[section]
%\newtheorem{lemma}{Lemma}[section]
%\newtheorem{proposition}{Proposition}[section]
%\newtheorem{definition}{Definition}
%\newtheorem{remark}{Remark}
%\newtheorem{cor}{Corollary}
\def\bep{\bs{\epsilon}}
\def\bL{\mb{L}}
\def\bsa{\mb{s}_a}
\def\bsb{\mb{s}_b}
\def\bsc{\mb{s}_c}
\def\bsi{\mb{s}_1}
\def\bso{\mb{s}_0}
\def\bY{\mb{Y}}
\def\bdelta{\bs{\delta}}
\def\balpha{\bs{\alpha}}
\def\bbeta{\bs{\beta}}
\def\bmu{\bs{\mu}}
\def\trans{^{\mbox{\sf \tiny T}}}
%\pagestyle{plain}
%\doublespacing
\title{Mendelian randomization using semiparametric linear transformation models}
\author{Yen-Tsung Huang}
\authormark{Huang}
\address{\orgdiv{Institute of Statistical Science}, \orgname{Academia Sinica}, \orgaddress{\state{Taipei}, \country{Taiwan}}}
\corres{Yen-Tsung Huang, No. 128 Academia Road Section 2, Nankang, Taipei 11529, Taiwan. \email{ythuang@stat.sinica.edu.tw}}
\presentaddress{}
\begin{document}
\abstract[Summary]{
Mendelian randomization (MR) uses genetic information as an instrumental variable (IV) to estimate the causal effect of an exposure of interest on an outcome in the presence of unknown confounding. We are interested in the causal effect of cigarette smoking on lung cancer survival, which is subject to confounding by underlying pulmonary functions. Despite the well-developed IV analyses for the continuous and binary outcomes, the scarcity of methodology for the survival outcome limits its utility for the time-to-event data collected in many observational studies. We propose an IV analysis method in the survival context, estimating causal effects on a transformed survival time and survival probabilities using semiparametric linear transformation models. We study the conditions under which hazard ratio and the effect on survival probability can be approximated. For statistical inference, we construct estimating equations to circumvent the difficulty in deriving joint likelihood of the exposure and the outcome, due to the unknown confounding. Asymptotic properties of the proposed estimators are established without parametric assumptions about confounders. We study the finite sample performance in extensive simulation studies. The MR analysis of a lung cancer study suggests a harmful prognostic effect of smoking packyears that would have been missed by the crude association.}
\keywords{estimating equations, hazard ratio, instrumental variable, semiparametric probit model, survival analysis}
\maketitle
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Introduction}
\label{ss:intro}
Randomized controlled trials (RCTs) are the gold standard for assessing the causal effect of an intervention on an outcome. Nevertheless, exposures like smoking ($Z$) are not suitable for RCTs, which renders observational studies a practical alternative. The bias due to unmeasured confounding in non-randomized observational studies can be of concern. Mendelian randomization (MR) provides a potential remedy to the confounding issue in non-randomized studies. \cite{katan, davey_smith} MR is an application of the instrumental variable (IV) method. \cite{didelez_sheehan}
MR uses genetic information ($G$) as an IV, which is assumed to be randomly assigned to the study subjects. MR analyses can be viewed as analyzing two RCTs: one for $G$-$Z$ (IV-exposure) association via a model of the exposure given the IV, $[Z|G]$, and the other one for $G$-$T$ (IV-outcome) via a model of survival time given the IV, $[T|G]$, to recover the causal effect of $Z$ on $T$ in the unknown model of survival time given $U$ and $Z$, $[T|U,Z]$ where $U$ indicates unknown confounders. IV analyses are well-developed for continuous and binary outcomes under regression models with certain link functions. \cite{lawlor, palmer} Recently, MR has also become a popular application in epidemiologic studies. \cite{granell_mr, ahmad_mr}
This paper is motivated by a lung cancer study where the role of cigarette smoking ($Z$) in lung cancer prognosis, i.e., time to death ($T$) is investigated. Despite the solid evidence of smoking in the etiology of lung cancer, \cite{smoking} the prognostic effect of smoking on lung cancer survival remains less well-established. Ferketich et al. (2013) \cite{ferketich} ($n$=4,200) and Bhatt et al. (2015) \cite{bhatt} ($n$=61,440) reported in two large studies that patients who were never smokers at the time of diagnosis had better survival compared with current smokers. However, given the fact that 80\% to 90\% of death in lung cancer patients is contributed by smoking, \cite{smoking} it is of importance to interrogate whether lung cancer prognosis is affected by intensity of smoking such as smoking pack-years. Smoking intensity may be influenced by socio-demographic factors, physiological conditions or underlying pulmonary functions ($U$) prior to the diagnosis of lung cancer, which are also determinants of patients' mortality. The detail regarding pre-clinical health conditions is not available in these studies. Even with adjustment of other observable factors such as age, gender and medical history, the association of smoking at diagnosis with lung cancer survival in standard analyses may still be subject to the undue influence by the underlying health condition, i.e., the so-called confounding. How to estimate the causal effect in a survival context is the goal of this study.
IV methodology for survival outcomes with right censoring is not as widely available as that for continuous or binary outcomes. In the context of clinical trials, Loeys and Goetghebeur (2003) \cite{loeys} and Yu et al. (2015) \cite{yu_jrssb} proposed correcting noncompliance for estimating true causal effects in a subset population, which can be construed as IV analyses limited to binary IV and exposure and may not be readily applicable to our problem with a continuous IV and exposure. MacKenzie et al. (2014) \cite{mackenzie} proposed an IV analysis for the hazard ratio (HR) using an additive-multiplicative hazard model. Li et al. (2015) \cite{li_iv} developed an IV method using Lin and Ying's additive hazard model \cite{lin_ying} to characterize the causal effect on the hazard difference; Tchetgen Tchetgen et al. (2015) \cite{tchetgen} and Martinussen et al. (2017) \cite{martinussen_2017} developed their methods using Aalen's additive hazard models \cite{aalen} to incorporate time-varying effects. Kjaersgaard and Parner (2016) \cite{kjaersgaard} proposed a pseudo-observation approach to IV analyses of cumulative incidence and restricted mean survival time. These existing methods focus on restricted model formulation, e.g., additive hazard model, or strong parametric assumptions due to the difficulty in characterizing corresponding effects in $[T|U,Z]$ by $[T|G]$ under many survival models.
Here we propose an IV analysis under linear transformation models, \cite{cheng_wei_ying} a broad class of models containing the widely used Cox proportional hazards (PH) models \cite{cox} and proportional odds (PO) models. \cite{bennett} The models relate a transformed survival time linearly to the exposure and confounders. The linearity well preserves collapsibility and thus proper causal interpretation. Under linear transformation models, our proposed method is closely related to the two stage IV methods and thus linked to the existing IV methodology for future generalization. Additionally, the model allows flexible parametric assumptions on the transformed survival time, which enables us to explore potential conditions where one may approximate various estimands of interest. Specifically, we introduce two causal effects: one on the expectation of the transformed survival time (Section~\ref{ss:effect_delta}), and the other on the entire survival probability (Section~\ref{ss:approx_omega}). For the effect on the expectation of the transformed survival time, we further study the condition under which the effect also approximates hazard ratio (HR) (Section~\ref{ss:approx_hr}). For statistical inference about our estimators, likelihood based approach is challenging due to the unknown confounding and thus we resort to a nonparametric estimating equation approach (Section~\ref{ss:infer}).
\section{The model and causal assumptions}
\label{ss:model}
Here we present two sets of models: data generating models in (\ref{tmodel0}) and (\ref{zmodel0}), and analysis models in (\ref{tmodel}) and (\ref{zmodel}) as well as required assumptions for IV analyses. In both data generating and analysis models, we introduce linear transformation models for the survival time and will later demonstrate its advantage over the hazard-based approach.
The goal of our proposed MR analysis is to investigate the causal effect of a continuous exposure $Z$ on the survival time $T$ by using a genetic determinant $G$ as an instrument. The association of smoking and the survival outcome is subject to confounding by an unknown factor $U$. The exposure is determined by a model of $Z$ given $G$ and $U$, $[Z|G,U]$, and the outcome is characterized by a survival model given $Z$ and $U$, $[T|Z,U]$. Because of censoring, $T$ is only observable up to a bivariate vector $(T^*, \delta)$, where $T^*=\min(T, C)$, $\delta=I(T\leq C)$, $I(\cdot)$ is an indicator function, and $C$ is the censoring time and is independent of $T$ given $Z$ and $G$.
The complete data consist of $n$ independent and identically distributed (iid) random vectors $\{T^*_i, \delta_i, Z_i, G_i, U_i\}$, where $i$ indexes subjects. We assume the following two data generating models for $[T|Z,U]$ and $[Z|G,U]$, respectively,
\begin{align}
H(T_i)&=-(\beta_0+\beta_ZZ_i+\beta_UU_i)+\epsilon_{Ti}\label{tmodel0}\\
Z_i&=\alpha_0+\alpha_GG_i+\alpha_UU_i+\epsilon_{Zi}\label{zmodel0},
\end{align}
where $U_i=\mu_U+\epsilon_{Ui}$; $\epsilon_{Ui}$, $\epsilon_{Ti}$ and $\epsilon_{Zi}$ follow arbitrary distributions with mean $0$ and respective variances $\sigma_U^2$, $\sigma_T^2$ and $\sigma_Z^2$ that are independent of other parameters and each other; $H(\cdot)$ is an unspecified strictly increasing smooth transformation function; $(\beta_0, \beta_U, \beta_Z)$ are unknown regression parameters representing the intercept, the effects of $U$ and $Z$ on the survival time $T$; and $(\alpha_0, \alpha_U, \alpha_G)$ are the regression parameters for the intercept, the effects of $U$ and $G$ on the exposure $Z$. Our parameter of interest is $\beta_Z$.
Several widely used survival models fall into the broad class of the semiparametric transformation model in (\ref{tmodel0}), with different parametric assumptions for $\epsilon_T$. Cox PH models \cite{cox} correspond to model (\ref{tmodel0}) with $\epsilon_T$ following the extreme value (EV) distribution, $P(\epsilon_Tt)$ by the above two models.
We characterize the causal effect and study identifiability assumptions using counterfactuals. \cite{rubin_78} Let $T(z)$ be the counterfactual survival time $T$ that would have been observed had $Z$ been set to $z$. $A\indep B|C$ denotes the independence of $A$ and $B$ conditional on $C$, and $A\notindep B|C$ denotes their conditional dependence. The following assumptions are sufficient conditions for Mendelian randomization analyses:
\begin{enumerate}
\item $G \notindep Z$: $G$ is associated with $Z$. A strong association ensures a strong IV. \label{a0}
\item $G\indep U$: the IV, $G$, is independent of the unknown confounder $U$. \label{a1}
\item $T\indep G|(Z, U)$: the survival time $T$ is independent of the IV ($G$) conditional on the exposure $Z$ and the unknown $U$. This is sometimes referred to as the exclusion restriction \cite{vanderweele_mr} or no pleiotropic effect of $G$ on $T$. \cite{lawlor} \label{a2}
\item $T(z)\indep Z|U$: there is no other confounding for the association of the survival time $T$ and the exposure $Z$ conditional on $U$. \label{a3}
\end{enumerate}
The assumptions \ref{a1} and \ref{a3} imply that there is no confounding for the association of the survival time $T$ and the instrument $G$, i.e., $T(g)\indep G$. By assumption~\ref{a1}, we have
\[
a_G=\alpha_G,
\]
$a_0=\alpha_0+\alpha_U\mu_U$, and $e_{Zi}=\alpha_U\epsilon_{Ui}+\epsilon_{Zi}$, where $e_Z$ is independent of $G$. $a_0$ and $a_G$ can be easily estimated by the least square estimators, $\hat{a}_0$ and $\hat{a}_G$.
To incorporate known confounding, the assumptions~\ref{a1}-\ref{a3} can be generalized to independence conditional on the measured confounders $\bW$: $G\indep U |\bW$, $T\indep G|(Z, U, \bW)$ and $T(z)\indep Z|(U, \bW)$. Suppose that confounding is linearly in $\bW$ and is not affected by $Z$, $G$ or $U$, models (\ref{tmodel}) and (\ref{zmodel}) can be modified accordingly by replacing $b_0$ and $a_0$ with $\bb_W\trans\bW$ and $\ba_W\trans\bW$, respectively, where the first element of $\bW$ is 1 for the intercept.
\section{Causal effects on a transformed survival time $H(T)$}
\label{ss:effect_delta}
Here we start with defining the average causal effect on $H(T)$ denoted by $\Delta$. By exploiting the linearity in (\ref{tmodel0}) and (\ref{tmodel}), we derive the expression of $\Delta$ by the parameters in data analysis models (\ref{tmodel}) and (\ref{zmodel}). We further show that its expression mimics the conventional two-stage IV method.
The average causal effect of $Z$ on $H(T)$ comparing $Z=z_1$ vs. $z_0$ is defined as $\Delta= \E[H(T(z_1))-H(T(z_0))]$, where $H(T(z))$ is a transformation $H(\cdot)$ of a counterfactual outcome $T$ had all subjects been assigned to receive the exposure $Z=z$. We show that
\begin{align*}
\E[H(T(z))|U]=\int H(t) dF_{T(z)}(t|U)
=\int H(t) dF_{T(z)}(t|U, z)=\int H(t) dF_T(t|U, z).
\end{align*}
The second equality is by assumption~\ref{a3} and the third equality is by consistency. \cite{vanderweele_consistency} With models (\ref{tmodel0}) and (\ref{zmodel0}),
\begin{align}
\Delta=\E\left\{\E[H(T(z_1))|U]\right\}-\E\left\{\E[H(T( z_0))|U]\right\}=-\beta_Z(z_1-z_0)\label{delta1}.
\end{align}
We introduce another counterfactual outcome $H(T(g))$, the outcome $H(T)$ had all subjects been assigned to receive the value of the instrument $G=g$. One can show that
\[
\E[H(T(g))]=\int H(t)dF_{H(T(g))}(H(t))=\int H(t)dF_{T(g)}(t)=\int H(t)dF_{T}(t|g)=\E[H(T)|G=g].
\]
The second last equality holds by $T(g)\indep G$ and consistency. It follows that
\begin{align}
\E[H(T(g))]&=\int H(t)dF_{T}(t|g)\nonumber\\
&=\int H(t) d\int\int F_T(t|Z, U)dF_Z(z|U, g)dF_U(u)\nonumber\\
&=-\left\{\beta_0+\beta_U\mu_U+\beta_Z(\alpha_0+\alpha_U\mu_U)+\beta_Z\alpha_Gg\right\}+\xi, \label{EHT}
\end{align}
where $\xi=\int\int\left\{-(\beta_U+\beta_Z\alpha_U)\epsilon_{Ui}-\beta_Z\epsilon_{Zi}\right\}dF_{\epsilon_Z}dF_{\epsilon_U}$ and if $\epsilon_Z$ and $\epsilon_U$ have mean zero, $\xi=0$. The second equality is by assumptions~\ref{a1} and \ref{a2}. As shown above, (\ref{EHT}) equals $\E[H(T)|G=g]$. We thus equate (\ref{EHT}) with the nonrandom part of model (\ref{tmodel}) and obtain $b_0=\beta_0+\beta_U\mu_U+\beta_Z(\alpha_0+\alpha_U\mu_U)+\xi$ and $b_G=\beta_Z\alpha_G$. Therefore,
\begin{align}
\Delta=-\frac{b_G}{a_G}(z_1-z_0) \label{delta2},
\end{align}
where $a_G\neq 0$ by assumption~\ref{a0}. Note that the results are general without any distributional assumption for $\epsilon_U$, $\epsilon_Z$ or $\epsilon_T$. Using counterfactuals, we show that the IV causal effect can be derived under the same assumptions (\ref{a0}-\ref{a3}) as the conventional IV method. We show in Section~\ref{ss:approx} that the effect in the scale of the hazard ratio or a difference in survival probabilities can be approximated under a normality assumption for $\epsilon_U$ and $\epsilon_Z$.
We denote $\E[H(T)|G]$ by $Y$ and rewrite as follows (assuming $\xi=0$):
\begin{align*}
Y&=-\left\{b_0^*+\beta_Z(a_0+a_GG)\right\},
\end{align*}
where $b_0^*=\beta_0+\beta_U\mu_U$. By plugging in $\hat{a}_0$ and $\hat{a}_G$, we have
\[
\hat{Y}=-\left\{b_0^*+\beta_Z\hat{Z}\right\},
\]
where $\hat{Z}=\hat{a}_0+\hat{a}_GG$. Because $\E[H(T)|G]=\E[\hat{Y}|G]$, $\beta_Z$ can be estimated by regressing $H(T)$ on $\hat{Z}$. Therefore, one can estimate $\beta_Z$ by a two-stage procedure similar to the conventional IV analyses: regress $Z$ on $G$ to obtain $\hat{Z}$ and then to regress $H(T)$ on $\hat{Z}$ to estimate the causal effect $\beta_Z$. $\Delta$ ($=-\beta_Z(z_1-z_0)$) in (\ref{delta2}) is equivalent to the two-stage procedure for $\beta_Z$. It follows that the residual inclusion procedure \cite{terza_2008} and the control function method \cite{wooldridge_2015} in the two stage estimation framework may also be incorporated in our method. Different from the commonly used two-stage least square estimator, $H(T)$ is not fully observable due to censoring. Statistical inference to be presented in Section~\ref{ss:infer} is developed in the context of survival outcome with censoring. Moreover, the above two-stage approach can only estimate the effect on $H(T)$. Additional development is required to characterize the effect on survival probability.
\section{Hazard ratio and survival probability}
\label{ss:approx}
In this section, we further demonstrate the advantage of linear transformation models of survival time. In addition to the effect on $\E[H(T)]$, flexibility of the model makes it possible to approximate hazard ratio and survival probability under certain conditions. In Section~\ref{ss:approx_hr}, we show that $e^{-\Delta}$ may approximate hazard ratio if we choose Cox PH or probit models in (\ref{tmodel}) along with other respective conditions. In Section~\ref{ss:approx_omega}, we study the conditions, under which the causal effect on the entire survival probability may be approximated and derive its closed form expression using the parameters in analysis models (\ref{tmodel}) and (\ref{zmodel}).
To exploit the conjugacy property of normal random variables, we make the following assumption,
\begin{enumerate}\setcounter{enumi}{4}
\item $\epsilon_U\sim N(0, \sigma_U^2)$ and $\epsilon_Z\sim N(0, \sigma_Z^2)$. \label{a4}
\end{enumerate}
With this additional parametric assumption, we show that under a Cox PH model in (\ref{tmodel0}), the hazard ratio can be approximated using a probit or a CoxPH model in (\ref{tmodel}) as the analysis model, and under a semiparametric probit model in (\ref{tmodel0}), an entire survival probability function can also be approximated using a probit model in (\ref{tmodel}) for analysis.
\subsection{Hazard ratio}
\label{ss:approx_hr}
Assuming that the data generating model (\ref{tmodel0}) is a Cox PH model and both the exposure and confounders are normally distributed, the error term $e_T$ in the resulting analysis model (\ref{tmodel}) is a convolution of EV and normal distributions. For estimating $\beta_Z$ (i.e., log hazard ratio), one has to derive the distribution of $e_T$ to specify the corresponding analysis model (\ref{tmodel}).
Under~\ref{a4}, $e_T$ in analysis model (\ref{tmodel}) equals $\epsilon_T+v$, a convolution of EV and normal distributions (EV-normal) with a probability density function (pdf)
\begin{align}
f_{e_T}(y)=\int_{-\infty}^{\infty}\exp\left(y-v-e^{y-v}\right)\frac{1}{\sqrt{2\pi}\sigma_v}\exp{\left(-\frac{v^2}{2\sigma_v^2}\right)}dv,
\label{ev_normal}
\end{align}
where $v=-(\beta_U+\beta_Z\alpha_U)\epsilon_{U}-\beta_Z\epsilon_{Z}\sim N(0, \sigma_v^2)$, $v$ and $\epsilon_T$ are independent by assumption~\ref{a3}, and $\sigma_v^2=(\beta_U+\alpha_U\beta_Z)^2\sigma_U^2+\beta_Z^2\sigma_Z^2$. We show in Figure~\ref{ff:ev_normal} that the EV-normal is approximated by a normal distribution under a large $\sigma_v^2$ (the blue lines), and that the EV-normal may only be approximated by an EV under a small $\sigma_v^2$ (i.e., small $\beta_Z$, $\alpha_U$ and $\beta_U$, the red lines). The simulation can be justified by the following convolution result as $\sigma_v^2\rightarrow\infty$. We express $f_{\epsilon_T}(y)$ as:
\begin{align*}
\int_{-\infty}^{\infty}\frac{1}{\sqrt{2\pi}\sigma_v}\exp{\left(-\frac{(y-u)^2}{2\sigma_v^2}\right)}\exp\left(u-e^{u}\right)du
&= \int_{-\infty}^{\infty}\frac{1}{\sqrt{2\pi}\sigma_v}\left(1-\frac{(y-u)^2}{2\sigma_v^2}\right)\exp\left(u-e^{u}\right)du+\rho_y\\
&=\frac{1}{\sqrt{2\pi}\sigma_v}\left(1-\frac{y^2+2\gamma y+\gamma^2+\pi^2/6}{2\sigma_v^2}\right)+\rho_y
\rightarrow \frac{1}{\sqrt{2\pi}\sigma_v}\exp\left(-\frac{(y+\gamma)^2}{2\sigma_v^2}\right),
\end{align*}
where $0<\rho_y<\bar{\rho}_y$, $\bar{\rho}_y=\int_{-\infty}^{\infty}\frac{1}{2\sqrt{2\pi}\sigma_v}\left(-\frac{(y-u)^2}{2\sigma_v^2}\right)^2\exp(u-e^u)du\rightarrow 0$ as $\sigma_v^2\rightarrow \infty$ because of the finite moments of $u$, and $\gamma$ is the Euler constant. A formal proof can be found in Cho and Huang (2019) \cite{cho_huang}. The above result shows that as $\sigma_v^2$ gets large, $e_T$ converges to a normal random variable with mean $-\gamma$ and variance $\sigma_v^2$. Practically, we found that the normal distribution with the variance of $e_T$, $\sigma_v^2+\pi^2/6$ has better approximation. The extra term $\pi^2/6$ arises from the exact moment matching and becomes negligible as $\sigma_v^2$ gets large. On the other hand, as $\sigma_v^2\rightarrow 0$, (\ref{ev_normal}) converges to $\exp(y-e^y)$, the pdf of EV.
Therefore, Figure~\ref{ff:ev_normal} and the above heuristic justification suggest that if the data generating model (\ref{tmodel0}) is a Cox PH model and both the exposure and confounders are normally distributed, the error term $e_T$ in the resulting analysis model (\ref{tmodel}) is a convolution of EV and normal distributions, which converges to a normal distribution when $\sigma_v^2$ gets large, or to a extreme value distribution when $\sigma_v^2$ is very small. The utility of this technical result is that fitting a probit model in (\ref{tmodel}) may approximate the Cox PH model in (\ref{tmodel0}) if the effect of interest or confounding is large; and if the effect and confounding are both very small, analysis with Cox PH model in (\ref{tmodel}) may approximate hazard ratio by $e^{-\Delta}$. We will further demonstrate this via MR analyses in Section~\ref{ss:sim}.
\subsection{Effects on survival probability $\mS_T(t)$}
\label{ss:approx_omega}
Assuming that the data generating model (\ref{tmodel0}) is a probit model (i.e., $\epsilon_T\sim N(0,1)$) and both the exposure and confounders are normally distributed (Assumption~\ref{a4}), the error term $e_T$ will also be normally distributed and thus the resulting analysis model (\ref{tmodel}) is another probit model. Under the probit model in (\ref{tmodel0}), the effect on the survival probability will be derived; and by fitting another probit model in (\ref{tmodel}), we show in Proposition~\ref{lem:reexp} the effect can be approximated.
Normality for $\epsilon_U$, $\epsilon_T$ and $\epsilon_Z$ appears to be a reasonable assumption for the lung cancer study (see Section~\ref{ss:app}). The resulting models (\ref{tmodel}) and (\ref{zmodel}) have
$e_{Zi}=\alpha_U\epsilon_{Ui}+\epsilon_{Zi}\sim N(0, \sigma_Z^{*2})$ and $e_{Ti}=-(\beta_U+\beta_Z\alpha_U)\epsilon_{Ui}-\beta_Z\epsilon_{Zi}+\epsilon_{Ti}\sim N(0, \sigma_T^{*2})$, where $\sigma_Z^{*2}=\alpha_U^2\sigma_U^2+\sigma_Z^2$ and $\sigma_T^{*2}=\sigma_T^2+\beta_Z^2\sigma_Z^2+(\beta_U+\beta_Z\alpha_U)^2\sigma_U^2$. We define the effect that compares the survival probability of $Z=z_1$ vs. $z_0$ using counterfactuals:
\begin{align}
\Psi(t|U)\equiv \mS_{T(z_1)}(t|U)-\mS_{T(z_0)}(t|U), \label{psi}
\end{align}
where $\mS_{T(z)}(t|U)=\mS_T(t|U,Z=z)$ by assumption \ref{a3} and consistency. A semiparametric probit model has the survival function:
\[
\mS_T(t|Z,U)=1-\Phi\left(\left[\log\Lambda(t)+\beta_0+\beta_UU+\beta_ZZ\right]\sigma_T^{-1}\right),
\]
where $\Lambda(t)=e^{H(t)}$ and $\Phi(\cdot)$ is the cumulative distribution function of a standard normal. The effect in (\ref{psi}) can be expressed as
\begin{align}
\Psi(t|U)=\mS_T(t|z_1,U)-\mS_T(t|z_0,U).
\end{align}
We prove the following result in Appendix.
\begin{proposition}
Suppose that $\epsilon_T\sim N(0, \sigma_T^2)$ and the normality assumption~\ref{a4} holds for $\epsilon_Z$ and $\epsilon_U$. If there exists a sequence indexed by $k=1,2,...$ of $\sigma_{Tk}^{*2}$ such that $\sigma_{Tk}^{*2}\rightarrow\sigma_T^2+\beta_Z^2\sigma_Z^2+(\beta_Z\alpha_U)^2\sigma_U^2$ as $k\to\infty$, then $\mS_T(t|U,z)\rightarrow 1-\Omega_z(t;\delta)$, where $\delta=\beta_U(U-\mu_U)$ and
\[
\Omega_z(t;\delta)=\Phi\left(\left[\log\Lambda(t)+b_0+\delta-\frac{b_G}{a_G}a_0+\frac{b_G}{a_G}z\right]\left(\sigma_{Tk}^{*2}-(\frac{b_G}{a_G})^2\sigma_Z^{*2}\right)^{-1/2}\right).
\]
\label{lem:reexp}
\end{proposition}
Denote $\Omega(t)=\Omega(t;\delta)\equiv\Omega_{z_0}(t; \delta)-\Omega_{z_1}(t;\delta)$ where $\Omega(t)$ is the limit that the effect on survival probability $\Psi(t|U)$ converges to under the assumptions stated in Proposition~\ref{lem:reexp}. The above result enables us to approximate $\Psi(t|U)$ by $\Omega(t;\delta)$ at various $U$ ($=\mu_U+\delta/\beta_U$) using the parameters in models (\ref{tmodel}) and (\ref{zmodel}):
\[\Psi(t|U)\approx -\left(\Omega_{z_1}(t;\delta)-\Omega_{z_0}(t;\delta)\right).\]
Under $\delta=0$, $\Omega(t;\delta=0)$ approximates the effect on survival probabilities conditional on the average value of $U$, i.e., $\Psi(t|\mu_U)=\mS_{T(z_1)}(t|\mu_U)-\mS_{T(z_0)}(t|\mu_U)$.
There are a few sufficient conditions for the condition $\sigma_{Tk}^{*2}\rightarrow\sigma_T^2+\beta_Z^2\sigma_Z^2+(\beta_Z\alpha_U)^2\sigma_U^2$ to hold. First, the confounder $U$ has a much larger effect on the outcome mediated through the exposure of interest $Z$, $\beta_Z\alpha_U$ than the effect independent of $Z$, $\beta_U$. Secondly, the true causal effect $\beta_Z$ is large. Thirdly, the variability of $U$ is small. Any of the above is sufficient to satisfy the condition. On the other hand, if the effect of interest is null, $\beta_Z=0$, the assumption is less likely to hold. For real applications, we suggest that the analysis of $\Psi(t)$ may be reserved as a secondary data analysis, and one may investigate $\Psi(t)$ only when the null $\Delta=0$ can be rejected and $\hat{\Delta}$ is reasonably large.
In summary, with normality of $\epsilon_Z$ and $\epsilon_U$, $-\log$ HR can be approximated by $\Delta$ when the underlying true model is a Cox PH model, and the survival probability can be approximated by $\Omega_z(t;\delta)$ when the survival time follows a semiparametric probit model. Notably, the normality assumption~\ref{a4} ensures additional interpretations but is not required for the validity of $\Delta$ as a causal effect on the transformed survival time or the following statistical inference.
\section{Inference about $\Delta$ and $\Omega(\cdot)$}
\label{ss:infer}
In the preceding two sections, closed form expressions of the causal effects $\Delta$ and $\Omega(t)$ have been derived with use of parameters in (\ref{tmodel}) and (\ref{zmodel}), the effect estimators $\hat{\Delta}$ and $\hat{\Omega}(t)$ can be easily constructed by plugging in estimators of the associated parameters: $\hat{\Delta}=\Delta(\hat{b}_G, \hat{a}_G)$ and $\hat{\Omega}(t)=\Omega(t; \hat{\Lambda}(t), \hat{b}_0, \hat{b}_G, \hat{a}_0, \hat{a}_G, \hat{\sigma}_Z^{*2})$. Here we estimate these parameters using a non-parametric estimating equation and study asymptotic normality of the two effect estimators $\hat{\Delta}$ and $\hat{\Omega}(t)$.
Due to the unknown confounding by $U$, we are not able to derive a joint likelihood of models (\ref{tmodel}) and (\ref{zmodel}). We instead construct estimating equations from the marginal likelihood of the two models. As will be shown, the inference imposes no particular distribution for $U$. Again, the parametric assumptions in Section~\ref{ss:approx} are for interpretations of $\Delta$ and $\Omega(t)$ as a $-\log$ HR and a difference in survival probabilities, respectively, but not for inference.
The observed data consist of $n$ iid random verctors $\{N_i(t), Y_i(t), \bXi\}$, where $N_i(t)=I(T^*\leq t, \delta=1)$ is the counting process for the observed failure event, $Y_i(t)=I(T_i^*\geq t)$ is the at-risk process, and $\bXi=(1, G_i)\trans$.
The log-likelihood function over a finite interval $t\in[0, \tau]$ \cite{zeng_lin_06} is:
\begin{align*}
\sum_{i=1}^n\left[\int_0^{\tau}\log\lambda(t)dN_i(t)
+\int_0^{\tau}\log\mathcal{G}'\left\{\int_0^tY_i(s)e^{b\trans\bXi}d\Lambda(s)\right\}dN_i(t)\right.\\
\left.+\int_0^{\tau}b\trans\bXi dN_i(t)-\mathcal{G}\left\{\int_0^{\tau}Y_i(t)e^{b\trans\bXi}d\Lambda(t)\right\}
\right],
\end{align*}
where $\lambda(t)=\frac{d\Lambda(t)}{dt}$; $\mathcal{G}'$ and $\mathcal{G}''$, respectively, are the first and second derivatives of $\mathcal{G}$. And, $\mathcal{G}(\cdot)$ is a known function, such as Box-Cox transformations, a class of logarithmic transformations, \cite{chen_jin_ying} a $-\log$ of survival function of a log normal random variable for probit models \cite{huang_cai} with mean zero and a known variance $\sigma_{T^*}^2$, and an identity function for Cox PH models. Similar to Zeng and Lin (2006), \cite{zeng_lin_06} we allow $\Lambda(\cdot)$ to be discrete, $\Lambda=(\Lambda_1, ...,\Lambda_m)\trans$, and replace $\lambda(t)$ with the jump size of $\Lambda$ at the event time $T_{(j)}^*$, $j=1,...,m$, $m$ is the number of events: $\lambda=(\lambda_1, ...,\lambda_m)\trans$, where $\{T_{(1)}^*, ..., T_{(m)}^*\}$ is a subset of $\{T_1^*, ..., T_n^*\}$ with $T_{(j)}^*=T_{(j)}$. We define the following notations: $a=(a_0, a_G)\trans$, $b=(b_0, b_G)\trans$, $\theta_a=(a\trans, \sigma_Z^{*2})\trans$, $\theta_b=(b\trans, \lambda\trans)\trans$, $\theta=(\theta_a\trans, \theta_b\trans)\trans$. A nonparametric maximum likelihood estimator (NPMLE) can be expressed as the following estimating equations for $\theta_b$, $U_{\theta_b}=(U_b\trans, U_{\Lambda}\trans)\trans$:
\begin{align*}
U_{b}=\sum_iU_{b i},\hspace{2em} U_{\Lambda_j}=\sum_iU_{\Lambda_j i},\hspace{5em} j=1,...,m,
\end{align*}
and
\begin{align*}
U_{b i}=&\left[\int_0^{\tau}dN_i(t)-\int_0^{\tau}\frac{\mathcal{G}''(t)}{\mathcal{G}'(t)}\int_0^tY_i(s)e^{b\trans\bXi}d\Lambda(s)dN_i(t)
+\mathcal{G}'(\tau)\int_0^{\tau}Y_i(t)e^{b\trans\bXi}d\Lambda(t)
\right]\bXi\\
U_{\Lambda_ji}=&\int_0^{T_{(j)}^*}dN_i(t)+\left\{\int_0^{\tau}\frac{\mathcal{G}''(t)}{\mathcal{G}'(t)}dN_i(t)-\mathcal{G}'(\tau)\right\}\int_0^{T_{(j)}^*}Y_i(t)e^{b\trans\bXi}d\Lambda(t).
\end{align*}
where $\mathcal{G}(t)=\mathcal{G}\left\{\int_0^tY_i(s)e^{b\trans\bXi}d\Lambda(s)
\right\}$. The least square estimator (LSE) for $\theta_a$ can be expressed as the following estimating equations denoted by $U_{\theta_a}=(U_a\trans, U_{\sigma_Z^{*2}})\trans$:
\begin{align*}
U_{a}=\sum_iU_{a i},\hspace{2em} U_{\sigma_Z^{*2}}=\sum_iU_{\sigma_Z^{*2} i},
\end{align*}
where $U_{a i}=(Z_i-\bXi\trans a)\bXi$ and $U_{\sigma_Z^{*2} i}=\sigma_Z^{*2}-(Z_i-\bXi\trans a)^2$.
Let $\hat{\theta}$ be the estimator solving the above estimating equations $U_{\theta}=(U_{\theta_a}\trans, U_{\theta_b}\trans)\trans=\mb{0}$. Since $\Delta$ and $\Omega(t)$, respectively are a function and a functional of $\theta$, we construct consistent estimators $\hat{\Delta}=\Delta(\hat{\theta})$ and $\hat{\Omega}(t)=\Omega(t; \hat{\theta})$ by plugging in $\hat{\theta}=(\hat{\theta}_a\trans, \hat{\theta}_b\trans)\trans$. The variance of $\hat{\Delta}$ and the covariance of $\hat{\Omega}(t)$ involve a joint distribution of $\hat{\theta}_a$ and $\hat{\theta}_b$, derived in Sections 1 and 2 of the Supplementary Materials. We present the following asymptotic results for $\Delta(\hat{\theta})$.
\begin{proposition}\label{thm:delta}
Suppose that $\theta_0$ is the true value of $\theta$, $\Delta$ is differentiable at $\theta$ and $U_{\theta}(\hat{\theta})=0$. As the sample size $n\rightarrow\infty$, $\sqrt{n}\left(\Delta(\hat{\theta})-\Delta(\theta_0)\right)$ converges to a normal distribution with mean 0 and variance
$\sigma_{\Delta}^2=D_\Delta\text{Cov}(U_{\theta i}(\theta_0))
D_{\Delta}\trans$,
where
$D_{\Delta}=\frac{\partial\Delta(\theta_0)}{\partial\theta}
\left[\text{E}\left(\frac{\partial U_{\theta i}(\theta_0)}{\partial \theta}\right)\right]^{-1}$.
\end{proposition}
For implementation, $\text{E}\left(\frac{\partial U_{\theta i}(\theta_0)}{\partial \theta}\right)$ can be estimated consistently by an empirical estimator $\frac{1}{n}\sum_{i=1}^{n}\frac{\partial U_{\theta i}(\hat{\theta})}{\partial \theta}=\frac{1}{n}\frac{\partial U_{\theta}(\hat{\theta})}{\partial \theta}$, and $\text{Cov}(U_{\theta i}(\theta_0))$ by the sample covariance of $U_{\theta i}(\hat{\theta})$. The derivatives of the estimating equations, $\frac{\partial U_{\theta}}{\partial \theta}$, is a block diagonal matrix of $\frac{\partial U_{\theta_a}}{\partial \theta_a}$ and $\frac{\partial U_{\theta_b}}{\partial \theta_b}$ with their respective expression provided in Section 4 of the Supplementary Materials.
\begin{proposition}\label{thm:omega}
Suppose that $\Omega(t)$ is compactly differentiable at $\theta$ and the normality assumptions in Proposition~\ref{lem:reexp} hold. As the sample size $n\rightarrow \infty$, $\sqrt{n}\left(\Omega(t; \hat{\theta})-\Omega(t; \theta_0)\right)$ converges to a zero-mean Gaussian process with a covariance matrix
$\Sigma=D_{\Omega}(t)
\text{Cov}(U_{\theta i}(\theta_0))
D_{\Omega}(t)\trans$,
where $D_{\Omega}(t)=\frac{\partial\Omega(t; \theta_0)}{\partial\theta}
\left[\text{E}\left(\frac{\partial U_{\theta i}(\theta_0)}{\partial \theta}\right)\right]^{-1}$.
\end{proposition}
Denote $\Omega_z(t)$ with the following expression:
\begin{align*}
\Omega_z(t)&=\mS_{\epsilon^*}\left(e^{\left(\log\Lambda(t)+b_0-\frac{b_G}{a_G}a_0+\frac{b_G}{a_G}z\right)\sigma^{*-1}}\right)\\
&=e^{-\mG\left( e^{\eta\sigma^{*-1}} \right)},
\end{align*}
where $\sigma^{*2}=\sigma_T^{*2}-\left(\frac{b_G}{a_G}\right)^2\sigma_Z^{*2}$, $\eta=\log\Lambda(t)+b_0-\frac{b_G}{a_G}a_0+\frac{b_G}{a_G}z$, $\mG(\cdot)=-\log\mS_{\epsilon^*}(\cdot)$; $\mS_{\epsilon^*}$ is the survival function of $\epsilon^*$, $\log\epsilon^*\sim N(0,1)$. We provide the expression of each element in $\frac{\partial\Omega(t; \theta_0)}{\partial\theta\trans}$ in Section 5 of the Supplementary Materials.
\section{Simulation}
\label{ss:sim}
Extensive numerical experiments were conducted to study the performance of the proposed methods under the finite sample. We generated samples of size $n=200$, a half of which had the instrumental variable $G=0$ and the other half had $G=1$. An unknown confounder $U$ was generated from a normal distribution with mean 0 and variance 4. The exposure of interest $Z$ was generated from $G$ and $U$ by the model: $Z_i=0+\alpha_U\times U_i+0.8\times G_i+\epsilon_{Zi}$, where $\epsilon_{Zi}$ follows the standard normal. The survival time $T$ was generated by the regression model: $\log T_i=-1.0+U_{i}\beta_U+Z_i\beta_Z+\epsilon_{Ti}$.
We generated data from various models ("true model in (\ref{tmodel0})" in Table~\ref{tt:cover}), including 1) probit models with $\epsilon_T$ following normal distribution with mean 0 and variance 4 (probit), 2) logistic models with $\epsilon_T$ following a non-normal distribution $\sqrt{12}logit(u)/\pi$ where $u\sim Uniform(0,1)$ (logistic), and 3) Cox PH models, i.e., with $\epsilon_T$ following $P(\epsilon_T95\%$ and lower power (Figure~\ref{ff:power}a). When the confounding effect was large ($\alpha_U=\beta_U=0.5$), the proposed CIs had a proper coverage if the effect of interest was not too large but the CP was slightly less than the expected under large effects. The performance was not sensitive to model misspecification when analyzing the survival time following logistic distribution using probit models (II in Table~\ref{tt:cover}), which shows the robustness of our proposed methods. Under the same causal effect, power decreased with the increase of confounding effects (Figures~\ref{ff:power}a and b).
We also conducted simulation studies where $\epsilon_T$ in (\ref{tmodel0}) was generated from the extreme value distribution, i.e., a Cox PH model and $\epsilon_Z$ was from the standard normal. We performed MR analyses with a Cox PH model in (\ref{tmodel}) even though the resulting $e_T$ followed a convolution of extreme value and normal distributions (EV-normal). Under such a misspecified model, we observed decreased CI coverages and large biases if the effect size $\beta_Z$ is large or confounding is strong (III in Table~\ref{tt:cover}). Probit models with $e_T\sim N(-\gamma, \sigma_v^2+\pi^2/6)$ performed well in terms of bias and the CP (IV in Table~\ref{tt:cover}) and had similar power, compared with Cox PH models (Figure~\ref{ff:power}d). $-\Delta$ obtained from MR with Cox PH models for (\ref{tmodel}) approximates the $\log$ HR under small effects and small confounding, and from probit models approximates the $\log$ HR under large effects or large confounding. The simulation studies confirm our theoretical results in Section~\ref{ss:approx_hr}.
The proposed estimator for $\Omega(t)$ and its confidence intervals covered the true $\Psi(t)$ well across various $(\alpha_U, \beta_U, \beta_Z)$ under a correctly specified model (Figure S1) and misspecified models (Figures S2 and S3). The results support the normal approximation for the EV-normal, shown in Figure~\ref{ff:ev_normal}.
For a simulation with 1000 replicates ($n$=200, 40\% censoring rate), the computation time using a desktop with the Intel\textsuperscript{\circledR} Core\textsuperscript{TM} i7-4790 CPU at 3.60 GHz/3.60 GHz and 16.00 GB RAM was 223.47 seconds and the convergence rate was 99.6\% on average.
\section{Smoking effect on lung cancer survival}
\label{ss:app}
A series of 148 snap-frozen tumor samples from stage I non-small cell lung cancer (NSCLC) patients with complete information on cigarette smoking was collected during surgery or biopsy from the Massachusetts General Hospital (MGH), Boston, MA, and the National Institute of Occupational Health, Oslo, Norway. \cite{huang_pnas} Demographic, smoking and other risk factor information was collected using a modified standardized American Thoracic Society respiratory questionnaire. \cite{zhou_ccr, huang_pnas} Among the 148 patients, the median and average follow up times were 3.9 and 4.5 years, respectively with censoring rate of 50\%. Genome-wide genotyping was carried out using Affymetrix 250K Nsp GeneChip (Affymetrix, Santa Clara, CA, USA).
Genotypes of each loci were coded as the additive mode (0/1/2), i.e., the number of minor allele, and the intensity of cigarette smoking was measured by the number of packs of cigarettes smoked per day multiplied by the number of years of smoking, i.e., pack-year, and square-root transformed on the collected pack-year value. For each genetic locus, we conducted multiple linear regression to assess its association with the square-root of smoking packyears, adjusting for age, sex, clinical stage (Ia vs. Ib), the population stratification and hospital. From a genome-wide analysis for the association with smoking packyears in 262,264 genetic loci (Figure~\ref{ff:manhattan}), we constructed a genetic score by the 20 loci: $G_i=\sum_{j=1}^{20}w_jS_{ij}$, where $S_{ij}$ is the number of the minor allele of genetic locus $j$ in subject $i$, and $w_j$ is the regression coefficient for the SNP-smoking association in the above linear regression. The genetic score captures considerable variability of cigarette smoking with $r^2$ of 0.49 (Figure S4). The 20 genetic loci was chosen with F-statistics greater than 10 for controlling the weak instrument bias under 1/10. \cite{burgess_thompson}
With the genetic score $G$, MR analyses were conducted to investigate the causal effect of cigarette smoking on lung cancer mortality using probit models and Cox PH models, adjusting for the covariates. Probit models showed that with a one unit increase in the square-root of smoking packyears, there was a significant decrease in the transformed survival time (Figure~\ref{ff:lcs}b), -0.33 (95\% CI=[-0.53, -0.12]; $p$-value=0.0019) (Table~\ref{tt:lcs}). The impact of smoking on the entire survival function conditional on the average value of $U$ is shown in Figure~\ref{ff:lcs}a, suggesting that smoking decreased the survival probabilities monotonically across the follow-up time. A one unit increase in the square-root of packyears from no smoking caused 3\% decrease in 5-year survival probability. Cox PH models showed that a one unit increase in smoking decreased $H(T)$ by 0.15 (95\% CI=[0.03, 0.27]; $p$-value=0.018) (Table~\ref{tt:lcs}). By assuming the magnitude of the effect and confounding is small, $\hat{\Delta}$ may approximate $-\log$ HR. Therefore, one may interpret the $e^{-\Delta}$ in MR Cox model as a hazard ratio of 1.16 (95\% CI=[1.03, 1.31]).
Compared to the conventional regression approach with and without adjustment of covariates (HR=1.14 and 1.06, respectively), MR analyses using probit and Cox PH models both reveal an even more harmful effect of cigarette smoking on lung cancer survival (Table~\ref{tt:lcs}). Similarly in Figure~\ref{ff:lcs}a, the effects on survival probabilities using the naive regression approach were smaller than those using MR analyses. The hazard difference (HD) estimated using additive hazard models \cite{tchetgen} showed a similar pattern (Table~\ref{tt:lcs}). The results support our hypothesis that deteriorating pulmonary functions prior to lung cancer diagnosis obscure the association of smoking at diagnosis with lung cancer survival. The confounding can be reduced but can not be completely eliminated by adjusting for collected covariates.
We have applied procedures and conducted sensitivity analyses to evaluate the exclusion restriction assumption~\ref{a2}. We excluded the genetic loci known to be associated with longevity, aging or pulmonary functions from the database that archived all published genetic loci for human traits in GWAS. \cite{gwas_catalog} We further conducted two sensitivity analyses. Mediation analyses for survival outcomes using probit models \cite{huang_cai} were conducted to investigate any effect of the genetic IV score on the survival not mediated through cigarette smoking (all $p>$0.2). In addition, we performed Egger regression \cite{egger} to examine the violation of the assumption with an intercept of 0.11 not significantly different from zero ($p$=0.33, Figure S5).
Model diagnostics were also performed to evaluate the goodness of fit of probit and Cox PH models. By comparing the largest departure of residuals of the probit models from the standard normal with their bootstrap counterparts, \cite{huang_cai} the diagnostics revealed no significant violation of the normality assumption ($p$-value=0.97). The PH assumption also held ($p$-value=0.47) by evaluating the weighted residuals. \cite{grambsch} The distribution of square root-transformed smoking packyears is approximately normal (Figure S6). Literatures support that the measured values of pulmonary function tests, i.e., the unknown confounder are normally distributed in a population. \cite{pft_book} The methods implemented in \textsc{Matlab} and data are available upon request.
\section{Discussion}
\label{ss:discuss}
In this article we propose an IV method in the survival context using flexible semiparametric linear transformation models. With standard assumptions~\ref{a0}-\ref{a3} for IV analyses, we derive the effect on a transformed survival time $\Delta$. With normality assumption in~\ref{a4}, we approximate the hazard ratio by $e^{-\Delta}$ using a Cox PH model or a probit model in (\ref{tmodel}), depending on the magnitude of the effect and confounding. If the true model (\ref{tmodel0}) is a probit model, we can further approximate the causal effect on survival probabilities. Although the normality assumptions provide additional interpretations, neither of validity and statistical inference about $\Delta$ depends on these assumptions.
The required assumptions are summarized in Table~\ref{tt:assum}. First, to interpret $\Delta$ as the effect on a transformed survival time $H(T)$ and to make statistical inference about $\hat{\Delta}$ and $\hat{\Omega}(t)$ only require Assumptions~\ref{a0}-~\ref{a3}. Secondly, to interpret $\Delta$ as $-\log$ HR has additional requirements: 1) Assumption~\ref{a4}, 2) data are generated from a Cox PH model, and 3a) under small treatment and confounding effects, a Cox PH model is used for data analysis, or 3b) under a large treatment or confounding effect, a probit model is used for data analysis. Thirdly, to interpret $\Omega(t)$ as the effect on survival probability additionally requires: 1) Assumption~\ref{a4}, 2) data are generated from a probit model and 3) a probit model is used for data analysis to estimate the associated parameters involved in $\Omega(t)$.
HR is a widely used effect estimate in medical applications, probably due to its proximity to risk ratio with relevant clinical utility and the popularity of Cox PH models. Existing methods may not estimate HR, \cite{martinussen_2017} or may approximate HR by very restricted models \cite{mackenzie} or by assuming a rare outcome. \cite{tchetgen, mackenzie} In our motivating example, the 5-year mortality rate of early stage lung cancer is not rare, ranging from 33 to 75\%. \cite{hoffman} The issue of causal interpretation of HR has also been discussed in literature. \cite{aalen_2015, hernan_hr} The hazard ratio we approximate corresponds to the causal hazard ratio defined by Aalen et al. (2015) \cite{aalen_2015}, which can be interpreted as controlled direct effect, the causal effect of the exposure on the survival at time $(t, t+\epsilon]$, conditional on unknown confounders $U$ and those surviving up to $t$; alternatively, it can also be construed as the causal effect of the exposure on the short time survival conditional on the unknown confounders $U$. Both provide valid causal interpretations; however, they do not coincide with the marginal hazard ratio in a randomized study due to non-collapsibility. The causal interpretation does not apply to the setting with a time-varying exposure. \cite{hernan_hr} The interpretation of the controlled direct effect emphasizes that the effect on the outcome is not through the unknown confounders $U$ since they have been controlled at a certain level and can no longer be influenced by the exposure. Whether the interpretation of HR as a controlled direct effect provide clinical or scientific relevance requires additional consideration.
%% Summary of assumptions
\section*{Supplementary Materials}
Derivations and proofs referenced in the text are provided in the eAppendix Section; results of
additional simulation (Figures S1-S3) and data analyses (Figures S4-S6) are also provided.
\section*{Appendix}
\paragraph{Proof of Proposition~\ref{lem:reexp}}
$\Psi(t)=\mS_{T}(t|\mu_U, Z=z_1)-\mS_{T}(t|\mu_U, Z=z_0)$, where $\mS_{T}(t|\mu_U, Z=z)=1-\Phi\left(\left[\log\Lambda(t)+\beta_0+\beta_U\mu_U+\beta_Zz\right]\sigma_T^{-1}\right)$. With the normality assumptions of $\epsilon_T$, $\epsilon_Z$ and $U$, model (\ref{tmodel}) can be shown to follow a semiparametric probit model:
\begin{align*}
\mS_{T}(t|G)&=\int\mS_T(t|G, U)dF_U(U|G)\\
&=\int\mS_T(t|G, U)dF_U(U)\q10uad (\text{by assumption}~\ref{a1})\\
&=\int\left(\int\mS_T(t|Z, U, G)dF_Z(z|U, G)\right)dF_U(U)\\
&=\int\left(\int\mS_T(t|Z, U)dF_Z(z|U, G)\right)dF_U(U)\qquad (\text{by assumption}~\ref{a2})\\
&=1-\Phi\left(\left[\log\Lambda(t)+\beta_0+\beta_U\mu_U+\beta_Z(\alpha_0+\alpha_GG+\alpha_U\mu_U)\right]\sigma_{Tk}^{*-1}\right)\\
&=1-\Phi\left(\left[\log\Lambda(t)+b_0+b_GG\right]\sigma_{Tk}^{*-1}\right).
\end{align*}
The second last equality is by assumption~\ref{a4}. It follows that
\begin{align*}
\Phi\left(\left[\log\Lambda(t)+\beta_0+\beta_UU+\beta_Zz\right]\sigma_T^{-1}\right)&=\Phi\left(\left[\log\Lambda(t)+b_0+\delta-\beta_Za_0+\beta_Zz\right]\sigma_T^{-1}\right)\\
&=\Phi\left(\left[\log\Lambda(t)+b_0+\delta-\frac{b_G}{a_G}a_0+\frac{b_G}{a_G}z\right]\sigma_T^{-1}\right),
\end{align*}
where $a_G\neq 0$ by assumption~\ref{a0}. Moreover, the converging value of $\sigma_{Tk}^{*2}$, denoted by $s_T^2$ can be expressed in the following,
\begin{align*}
s_T^{2}&=\sigma_T^2+\beta_Z^2\sigma_Z^2+(\beta_Z\alpha_U)^2\sigma_U^2\\
&=\sigma_T^2+(\frac{b_G}{a_G})^2\sigma_Z^{*2}.
\end{align*}
We denote $\eta_t=\log\Lambda(t)+\beta_0+\beta_UU+\beta_Zz$ in order to simplify the notation and re-express $\Phi\left(\left[\log\Lambda(t)+\beta_0+\beta_UU+\beta_Zz\right]\sigma_T^{-1}\right)$ by a Taylor expansion,
\begin{align*}
\Phi(\eta_t\sigma_T^{-1})&=\Phi\left(\eta_t\left\{s_T^2-\left(\frac{b_G}{a_G}\right)^2\sigma_Z^{*2}\right\}^{-1/2}\right)\\
&=\Phi\left(\eta_t\left\{\sigma_{Tk}^{*2}-\left(\frac{b_G}{a_G}\right)^2\sigma_Z^{*2}\right\}^{-1/2}\right)+\rho(s_T^2-\sigma_{Tk}^{*2}),
\end{align*}
where $\rho=-\frac{1}{2}\phi\left(\eta_t\left\{s_T^{*2}-\left(\frac{b_G}{a_G}\right)^2\sigma_Z^{*2}\right\}^{-1/2}\right)\left(s_T^{*2}-\left(\frac{b_G}{a_G}\right)^2\sigma_Z^{*2}\right)^{-2/3}$, $|s_T^{*2}-s_T^2|\leq|\sigma_T^{*2}-s_T^2|$ and $\phi(z)=\partial\Phi(z)/\partial z$.
We note that if the survival probability is not extremely low or high, then $\rho$ tends to be much smaller than 1. Consequently, $\rho(s_T^2-\sigma_{Tk}^{*2})\rightarrow 0$ faster than $s_T^2-\sigma_{Tk}^{*2}$. Therefore, we show that $\Phi\left(\left[\log\Lambda(t)+\beta_0+\beta_UU+\beta_Zz\right]\sigma_T^{-1}\right)$ with parameters in models (\ref{tmodel0}) and (\ref{zmodel0}) converges to $\Phi\left(\left[\log\Lambda(t)+b_0+\delta-\frac{b_G}{a_G}a_0+\frac{b_G}{a_G}z\right]\left(\sigma_T^{*2}-(\frac{b_G}{a_G})^2\sigma_Z^{*2}\right)^{-1/2}\right)$ with parameters in models (\ref{tmodel}) and (\ref{zmodel}).
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%
% References %
%%%%%%%%%%%%%%
%\bibliographystyle{WileyNJD-AMA}
\bibliography{MR_survival}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newpage
%% Simulation: Delta
\begin{table}[htp]
\caption{Simulation results of 1000 replicates. Coverage of 95\% CI and standard deviation, $SD$ are compared between the method assuming independence of $a_G$ and $b_Z$ ($V_{ind}$) and that accounting for their correlation (Adjust). }
\centering
\footnotesize
\begin{tabular}{c|ccccc|ccccc}\hline
\multirow{3}{*}{$\beta_Z$} & \multicolumn{5}{|c}{$\alpha_U=\beta_U=0.2$} & \multicolumn{5}{|c}{$\alpha_U=\beta_U=0.5$} \\\cline{2-11}
& Bias & \multicolumn{2}{c}{CI coverage (\%)} & \multicolumn{2}{c|}{$SD\times 100$} & Bias & \multicolumn{2}{c}{CI coverage (\%)} & \multicolumn{2}{c}{$SD\times 100$}\\
& $\times 100$ & Adjust & $V_{ind}$ & Adjust & $V_{ind}$ & $\times 100$ & Adjust & $V_{ind}$ & Adjust & $V_{ind}$\\\hline
& \multicolumn{10}{c}{I. Analysis model in (\ref{tmodel}): probit; true model in (\ref{tmodel0}): probit}\\\cline{2-11}
0.0 & 0.32 & 95.8 & 95.9 & 20.6 & 20.8 & -0.26 & 96.0 & 95.7 & 20.6 & 20.9 \\
0.2 & -0.29 & 95.2 & 95.5 & 20.5 & 21.0 & -2.21 & 95.7 & 96.6 & 20.3 & 21.1 \\
0.4 & -0.54 & 94.7 & 96.0 & 20.7 & 21.8 & -4.82 & 94.9 & 96.1 & 20.1 & 21.8 \\
0.6 & -1.04 & 95.2 & 96.9 & 21.0 & 23.0 & -7.94 & 92.8 & 95.7 & 20.1 & 23.0 \\
0.8 & -1.01 & 95.7 & 97.9 & 21.4 & 24.6 & -11.3 & 91.5 & 95.1 & 20.2 & 24.6 \\
1.0 & -0.51 & 95.1 & 98.2 & 22.0 & 26.5 & -14.6 & 88.6 & 95.0 & 20.6 & 26.4 \\
\hline
& \multicolumn{10}{c}{II. Analysis model in (\ref{tmodel}): probit; true model in (\ref{tmodel0}): logistic}\\\cline{2-11}
0.0 & -0.72 & 95.7 & 96.1 & 20.7 & 21.0 & -0.89 & 96.0 & 96.3 & 20.8 & 21.0 \\
0.2 & -0.07 & 95.7 & 96.2 & 20.7 & 21.2 & -1.95 & 95.8 & 96.6 & 20.4 & 21.2 \\
0.4 & 0.27 & 95.3 & 96.5 & 20.8 & 21.9 & -4.19 & 95.7 & 97.1 & 20.1 & 21.9 \\
0.6 & 0.72 & 96.0 & 97.6 & 21.1 & 23.2 & -7.05 & 95.8 & 97.0 & 20.1 & 23.1 \\
0.8 & 1.26 & 96.1 & 97.9 & 21.6 & 24.8 & -10.6 & 94.9 & 96.8 & 20.2 & 24.6 \\
1.0 & 1.88 & 96.0 & 98.6 & 22.1 & 26.7 & -14.0 & 92.1 & 96.6 & 20.6 & 26.5 \\
\hline
& \multicolumn{10}{c}{III. Analysis model in (\ref{tmodel}): Cox PH; true model in (\ref{tmodel0}): Cox PH }\\\cline{2-11}
0.0 & 0.14 & 96.8 & 97.0 & 29.6 & 29.7 & -2.61 & 97.1 & 98.1 & 31.7 & 31.4\\
0.2 & -2.07 & 97.7 & 97.2 & 29.1 & 29.9 & -7.21 & 98.0 & 98.5 & 29.9 & 31.3\\
0.4 & -5.08 & 96.8 & 97.3 & 28.3 & 30.5 & -15.4 & 97.6 & 96.9 & 28.3 & 31.6\\
0.6 & -10.1 & 95.4 & 96.5 & 27.3 & 31.4 & -25.4 & 92.4 & 93.1 & 26.9 & 32.2\\
0.8 & -18.0 & 90.7 & 94.3 & 26.4 & 32.3 & -37.5 & 70.1 & 84.1 & 25.7 & 32.8\\
1.0 & -28.0 & 78.7 & 88.2 & 25.5 & 33.2 & -51.2 & 39.0 & 67.4 & 24.6 & 33.4\\
\hline
& \multicolumn{10}{c}{IV. Analysis model in (\ref{tmodel}): probit; true model in (\ref{tmodel0}): Cox PH }\\\cline{2-11}
0.0 & 1.11 & 96.2 & 96.2 & 39.9 & 40.3 & -2.50 & 96.6 & 97.9 & 42.7 & 42.7\\
0.2 & 4.63 & 96.0 & 97.3 & 39.2 & 40.6 & -1.94 & 96.8 & 98.5 & 41.1 & 43.4\\
0.4 & 6.99 & 95.3 & 97.7 & 38.1 & 41.3 & -2.06 & 96.9 & 98.9 & 42.3 & 47.9\\
0.6 & 7.04 & 95.4 & 98.4 & 36.8 & 42.4 & -1.73 & 96.9 & 99.1 & 43.7 & 53.2\\
0.8 & 3.60 & 96.2 & 98.8 & 35.3 & 43.4 & -1.29 & 96.5 & 99.8 & 45.3 & 59.4\\
1.0 & 1.97 & 96.1 & 99.1 & 35.4 & 46.5 & -1.06 & 96.1 & 99.8 & 47.0 & 65.9\\
\hline
\end{tabular}
\label{tt:cover}
\end{table}
%% LCS: delta
\begin{table}[htp]
\caption{Effects on $H(T)$ by a one unit increase in square root of smoking packyears. MR, $V_{ind}$ indicates MR analyses assuming independence of $\hat{\theta}_a$ and $\hat{\theta}_b$; "Naive, adjusted" indicates standard analyses for the association of smoking and the outcome, adjusting for age, sex, clinical stage, cohort type and the first four principal components for population stratification as covariates; "Naive, unadjusted" indicates regression analyses without adjusting for covariates. $\Delta$ in Cox PH models can be interpreted as $-\log$HR}
\centering
\begin{tabular}{clccc}\hline
Model & Method & $\Delta$ or HD & 95\% CI & $p$-value \\\hline
\multirow{4}{2cm}{Probit} & MR & -0.33 & -0.53, -0.12 & 0.0019\\
& MR, $V_{ind}$ & -0.33 & -0.63, -0.03 & 0.034\\
& Naive, adjusted & -0.21 & -0.42, -0.01 & 0.042\\
& Naive, unadjusted & -0.10 & -0.26, 0.06 & 0.23\\\hline
\multirow{4}{2cm}{Cox PH} & MR & -0.15 & -0.27, -0.03 & 0.018\\
& MR, $V_{ind}$ & -0.15 & -0.30, -0.01 & 0.058\\
& Naive, adjusted & -0.13 & -0.24, -0.02 & 0.022\\
& Naive, unadjusted & -0.06 & -0.14, 0.03 & 0.19\\\hline
\multirow{3}{2cm}{Additive hazard} & MR & 3.68$\times 10^{-5}$ & 0.53$\times 10^{-5}$, 7.25$\times 10^{-5}$ & 0.020\\
& Naive, adjusted & 3.30$\times 10^{-5}$ & 0.84$\times 10^{-5}$, 5.77$\times 10^{-5}$ & 0.0086\\
& Naive, unadjusted & 1.82$\times 10^{-5}$ & -0.13$\times 10^{-5}$, 3.78$\times 10^{-5}$ & 0.068\\\hline
\end{tabular}
\label{tt:lcs}
\end{table}
\begin{table}[htp]
\caption{Required causal and model assumptions for interpretation and statistical inference of $\Delta$ and $\Omega(t)$.}
\centering
\begin{tabular}{r|ccc|cc}\hline
\multirow{3}{*}{Assumptions} & \multicolumn{3}{c|}{Interpretation} & \multicolumn{2}{c}{Inference}\\
& $\Delta=$ effect & $\Delta=$ & $\Omega(t)=$ effect & \multirow{2}{*}{$\Delta$} & \multirow{2}{*}{$\Omega(t)$}\\
& on $H(T)$ & $-\log$ HR & on $\mS_T(t)$ & & \\\hline
\ref{a0}-\ref{a3} & ${\surd}$ & ${\surd}$ & ${\surd}$ & ${\surd}$ & ${\surd}$\\
\ref{a4} & - & ${\surd}$ & ${\surd}$ & - & - \\
True model (\ref{tmodel0}) & - & Cox & Probit & - & - \\
Analysis model (\ref{tmodel}) under &&&&&\\
small $\beta_Z$, $\beta_U$, $\alpha_U$ & - & Cox & Probit & - & - \\
large $\beta_Z$, $\beta_U$, $\alpha_U$ & - & Probit & Probit & - & - \\\hline
\end{tabular}
\label{tt:assum}
\end{table}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newpage
%% EV-normal
\begin{figure}[htp]
\caption{Probability density of $e_T=\epsilon_T+v$, a convolution of an extreme value random variable $\epsilon_T$ and a normal random variable $v\sim N(0, \sigma_v^2)$ (solid lines). Dashed lines represent normal approximation and dotted lines represent the approximation by the extreme value distribution, $1-\exp(-e^{y})$, both by matching the first two moments. }
\centering
\includegraphics[angle=0, width=8cm, height=8cm]{ev_normal_pdf.eps}
\label{ff:ev_normal}
\end{figure}
%% Power: probit and logit
\begin{figure}[htp]
\caption{\% of $p<0.05$ under different configurations of true effects $\beta_Z$ and confounding effects $\alpha_U$ and $\beta_U$. $V_{ind}$ indicates analyses using the variance estimator assuming independence of $\hat{a}_G$ and $\hat{b}_G$ (dotted lines); black lines (solid and dotted) indicate no confounding ($\alpha_U=\beta_U=0$); red lines (solid and dotted) indicate weak confounding ($\alpha_U=\beta_U=0.2$); and blue lines (solid and dotted) indicate strong confounding ($\alpha_U=\beta_U=0.5$). }
\centering
\subfigure[]{\includegraphics[angle=0, width=7cm, height=7cm]{probit_power.eps}}
\subfigure[]{\includegraphics[angle=0, width=7cm, height=7cm]{logistic_power.eps}}
\subfigure[]{\includegraphics[angle=0, width=7cm, height=7cm]{Cox_power.eps}}
\subfigure[]{\includegraphics[angle=0, width=7cm, height=7cm]{probit_normapp_power.eps}}
\label{ff:power}
\end{figure}
%% GWAS of smoking
\begin{figure}[htp]
\caption{Genome-wide analyses of genetic association with smoking packyears. Red dots, the 20 genetic loci with $F$-statistic greater than 10. }
\centering
\includegraphics[angle=0, width=16cm, height=6cm, trim={0cm, 0cm, 0cm, 2cm}]{manhattan.eps}
\label{ff:manhattan}
\end{figure}
%% LCS: Omega
\begin{figure}[htp]
\caption{MR results of lung cancer study comparing $s=1$ vs. $0$. (a) Effect of smoking on lung cancer survival $\Omega(t)$ estimated by MR and regression (Naive analysis); (b) estimated transformation of survival time in probit models. $\Omega(t)$: the effect on survival probability comparing $s=1$ vs. $0$; $H(t)$: estimated nonparametric transformation of $t$. }
\centering
\subfigure[Smoking effects on survival probabilities]{\includegraphics[angle=0, width=7cm, height=7cm]{LCS_smoke_death.pdf}}
\subfigure[Estimated transformation function $H(t)$]{\includegraphics[angle=0, width=7cm, height=7cm]{LCS_tht.eps}}
\label{ff:lcs}
\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\end{document}
PKÒ}ŠíöíöPK
nTBWÐNÐN
MR_survival.bibUT ñœeñœe
@ARTICLE{pearl_dag,
author={Pearl, J.},
title={Causal diagrams for empirical research},
journal={Biometrika},
year={1995},
volume={82},
number={},
pages={669-709},
note={}
}
@ARTICLE{rubin_78,
author={Rubin, D.},
title={Bayesian inference for causal effects},
journal={Annals of Statistics},
year={1978},
volume={6},
number={},
pages={34-58},
note={}
}
@ARTICLE{lange,
author={Lange, T. and Hansen, J. V.},
title={Direct and indirect effects in a survival context},
journal={Epidemiology},
year={2011},
volume={22},
number={4},
pages={575-581},
note={}
}
@ARTICLE{huang_cai,
author={Huang, Y. T. and Cai, T.},
title={Mediation analysis for survival data using semiparametric probit models},
journal={Biometrics},
year={2016},
volume={72},
number={},
pages={563-574},
note={}
}
@ARTICLE{huang_epi17,
author={Huang, Y. T. and Yang, H. I.},
title={Causal mediation analysis of survival outcome with multiple mediators},
journal={Epidemiology},
year={2017},
volume={},
number={},
pages={},
note={in press}
}
%%%%%%%%%%%%%%%
@ARTICLE{ferketich,
author={Ferketich, A. K. and Niland, J. C. and Mamet, R. and Zornosa, C. and D'Amico, T. A. and Ettinger, D. S. and Kalemkerian, G. P. and Pisters K. M. and Reid, M. E. and Otterson, G. A.},
title={Smoking status and survival in the national comprehensive cancer network non-small cell lung cancer cohort},
journal={Cancer},
year={2013},
volume={119},
number={},
pages={847-53},
note={}
}
@ARTICLE{bhatt,
author={Bhatt, V. R. and Batra, R. and Silberstein, P. T. and Loberiza, F. R. Jr. and Ganti, A. K.},
title={Effect of smoking on survival from non-small cell lung cancer: a retrospective Veterans' Affairs Central Cancer Registry (VACCR) cohort analysis},
journal={Medical Oncology},
year={2015},
volume={32},
number={},
pages={339},
note={}
}
@ARTICLE{smoking,
author={{U.S. Department of Health and Human Services}},
title={The Health Consequences of Smoking: A Report of the {U.S. Surgeon General}},
journal={},
year={2004},
volume={},
number={},
pages={},
note={}
}
@ARTICLE{granell_mr,
author={Granell, R. and Henderson, A. J. and Evans, D. M. and Smith, G. D. and Ness, A. R. and Lewis, S. and Palmer, T. M. and Sterne, J. A. C.},
title={Effects of {BMI}, fat mass, and lean mass on asthma in childhood: a {Mendelian} randomization study},
journal={PLoS Medicine},
year={2014},
volume={11},
number={},
pages={e1001669},
note={}
}
@ARTICLE{ahmad_mr,
author={Ahmad, O. S. and Morris, J. A. and Mujammami, M. and Forgetta, V. and Leong, A. and Li, R. and Turgeon, M. and Greenwood, C. M. T. and Thanassoulis, G. and Meigs, J. B. and Sladek, R. and Richards, J. B.},
title={A {Mendelian} randomization study of the effect of type-2 diabetes on coronary heart disease},
journal={Nature Communications},
year={2015},
volume={6},
number={},
pages={7060},
note={}
}
@ARTICLE{zhang_mr,
author={Zhang, C. and Doherty, J. A. and Burgess, S. and Hung, R. J. and Lindstrom, S. and Kraft, P. and Gong, J. and Amos, C. I. and Sellers, T. A. Monteiro, A. N. A. and Chenevix-Trench, G. and Bickeboller, H. and Risch, A. and Brennan, P. and Mckay, J. D.
},
title={Genetic determinants of telomere length and risk of common cancers: a {Mendelian} randomization study},
journal={Human Molecular Genetics},
year={2015},
volume={24},
number={},
pages={5356-5366},
note={}
}
{and Houlston, R. S. and Landi, M. T. and Timofeeva, M. N. and Wang, Y. and Heinrich, J. and Kote-Jarai, Z. and Eeles, R. A. and Muir, K. and Wiklund, F. and GrÃ¶nberg, H. and Berndt, S. I. and Chanock, S. J. and Schumacher, F. and Haiman, C. A. and Henderson, B. E. and Olama, A. A. A. and Andrulis, I. L. and Hopper, J. L. and Chang-Claude, J. John, E. M. and Malone, K. E. and Gammon, M. D. and Ursin, G. and Whittemore, A. S. and Hunter, D. J. and Gruber, S. B. and Knight, J. A. and Hou, L. and Marchand, L. L. and Newcomb, P. A. and Hudson, T. J. and Chan, A. T. and Li, L. and Woods, M. O. and Ahsan, H. and Pierce, B. L.}
@ARTICLE{tchetgen,
author={Tchetgen Tchetgen, E. J. and Walter, S. and Vansteelandt, S. and Martinussen, T. and Glymour, M.},
title={Instrumental variable estimation in a survival context},
journal={Epidemiology},
year={2015},
volume={26},
number={},
pages={401-410},
note={}
}
@ARTICLE{aalen,
author={Aalen, O. O.},
title={A linear regression model for the analysis of life times},
journal={Statistics in Medicine},
year={1989},
volume={8},
number={},
pages={907-925},
note={}
}
@ARTICLE{doll,
author={Doll, R. and Hill, A. B.},
title={Smoking and carcinoma of the lung; preliminary report},
journal={British Medical Journal},
year={1950},
volume={2},
number={},
pages={739-948},
note={}
}
@ARTICLE{davey_smith,
author={Davey Smith, G. and Ebrahim, S.},
title={{'Mendelian randomisation': can genetic epidemiology contribute to understanding environmental determinants of disease?}},
journal={International Journal of Epidemiology},
year={2003},
volume={32},
number={},
pages={1-22},
note={}
}
@ARTICLE{katan,
author={Katan, M. B.},
title={{Apolipoprotein E isoforms, serum cholesterol, and cancer}},
journal={Lancet},
year={1986},
volume={1},
number={},
pages={507-508},
note={}
}
@book{bowden_turkington,
author = {Bowden, R. J. and Turkington, D. A.},
title = {Instrumental Variables},
publisher = {Cambridge University Press},
year = 1984,
address = {Cambridge, United Kingdom},
note = {},
isbn = {}
}
@ARTICLE{cox,
author={Cox, D. R.},
title={Regression models and life-tables. },
journal={Journal of the Royal Statistical Society, Series B},
year={1972},
volume={34},
number={},
pages={187-220},
note={}
}
@ARTICLE{andersen_gill,
author={Andersen, P. K. and Gill, R. D.},
title={{Cox's regression model for counting process: a large sample study}},
journal={Annals of Statistics},
year={1982},
volume={10},
number={},
pages={1100-1120},
note={}
}
@ARTICLE{bennett,
author={Bennett, S.},
title={Analysis of survival data by the proportional odds model},
journal={Statistics in Medicine},
year={1983},
volume={2},
number={},
pages={273-277},
note={}
}
@ARTICLE{cheng_wei_ying,
author={Cheng, S. C. and Wei, L. J. and Ying, Z.},
title={Analysis of transformation models with censored data},
journal={Biometrika},
year={1995},
volume={82},
number={},
pages={835-845},
note={}
}
@ARTICLE{chen_jin_ying,
author={Chen, K. and Jin, Z. and Ying, Z.},
title={Semiparametric analysis of transformation models with censored data},
journal={Biometrika},
year={2002},
volume={89},
number={},
pages={659-668},
note={}
}
@ARTICLE{zeng_lin,
author={Zeng, D. and Lin, D. Y.},
title={Maximum likelihood estimation in semiparametric regression models with censored data},
journal={Journal of the Royal Statistical Society, Series B},
year={2007},
volume={69},
number={},
pages={507-564},
note={}
}
@ARTICLE{zeng_lin_06,
author={Zeng, D. and Lin, D. Y.},
title={Efficient estimation of semiparametric transformation models for counting process},
journal={Biometrika},
year={2006},
volume={93},
number={},
pages={627-640},
note={}
}
@ARTICLE{hoffman,
author={Hoffman, P. C. and Mauer, A. M. and Vokes, E. E.},
title={Lung cancer},
journal={Lancet},
year={2000},
volume={355},
number={},
pages={479-485},
note={}
}
@ARTICLE{lawlor,
author={Lawlor, D. A. and Harbord, R. M. and Sterne, J. A. C. and Timpson, N. and Davey Smith, G.},
title={Mendelian randomization: using genes as instruments for making causal inference in epidemiology},
journal={Statistics in Medicine},
year={2008},
volume={27},
number={},
pages={1133-1163},
note={}
}
@ARTICLE{palmer,
author={Palmer, T. M. and Thompson, J. R. and Tobin, M. D. and Sheehan, N. A. and Burton, P. R.},
title={Adjusting for bias and unmeasured confounding in {Mendelian} randomization studies with binary responses},
journal={International Journal of Epidemiology},
year={2008},
volume={37},
number={},
pages={1161-1168},
note={}
}
@ARTICLE{didelez_sheehan,
author={Didelez, V. and Sheehan, N.},
title={{Mendelian} randomization as an instrumental variable approach to causal inference},
journal={Statistical Methods in Medical Research},
year={2007},
volume={16},
number={},
pages={309-330},
note={}
}
@ARTICLE{vanderweele_mr,
author={VanderWeele, T. J. and Tchetgen Tchetgen, E. J. and Cornelis, M. and Kraft, P.},
title={Methodological challenges in {Mendelian} randomization},
journal={Epidemiology},
year={2014},
volume={25},
number={},
pages={427-435},
note={}
}
@ARTICLE{burgess,
author={Burgess, S.},
title={Sample size and power calculations in {Mendelian randomization} with a single instrumental variable and a binary outcome},
journal={International Journal of Epidemiology},
year={2014},
volume={43},
number={},
pages={922-929},
note={}
}
@ARTICLE{lin_wang_probit,
author={Lin, X. and Wang, L.},
title={A semiparametric probit model for case 2 interval-censored failure time data},
journal={Statistics in Medicine},
year={2010},
volume={29},
number={},
pages={972-981},
note={}
}
@ARTICLE{lin_probit,
author={Lin, H. and Zhou, L. and Li, C. and Li, Y.},
title={Semiparametric transformation models for semicompeting survival data},
journal={Biometrics},
year={2014},
volume={70},
number={},
pages={599-607},
note={}
}
@techreport{iarc_x,
title = {{World Cancer Report 2014}},
author = {{International Agency for Research on Cancer}},
group = {csg},
year = {2014},
institution = {{International Agency for Research on Cancer}},
month = {},
Date-Added = {},
Date-Modified = {}
}
@book{iarc,
author = {{IARC}},
title = {{World Cancer Report 2014}},
publisher = {World Health Organization},
year = 2014,
volume = {},
series = {},
address = {Lyon, France},
note = {},
isbn = {}
}
@ARTICLE{huang_pnas,
author={Huang, Y. T. and Lin, X. and Liu, Y. and Chirieac, L. R. and McGovern, R. and Wain, J. and Heist, R. and Skaug, V. and Zienolddiny, S. and Haugen, A. and Su, L. and Fox, E. A. and Wong, K. K. and Christiani, D. C.},
title={Cigarette smoking increases copy number alterations in nonsmall-cell lung cancer},
journal={Proceedings of the National Academy of Sciences},
year={2011},
volume={108},
number={},
pages={16345-16350},
note={}
}
@ARTICLE{grambsch,
author={Grambsch, P. and Therneau, T.},
title={Proportional hazards tests and diagnostics based on weighted residuals},
journal={Biometrika},
year={1994},
volume={81},
number={},
pages={515-526},
note={}
}
@ARTICLE{vanderweele_consistency,
author={VanderWeele, T. J.},
title={Concerning the consistency assumption in causal inference},
journal={Epidemiology},
year={2009},
volume={20},
number={},
pages={880-883},
note={}
}
@ARTICLE{li_iv,
author={Li, J. and Fine, J. and Brookhart, A.},
title={Instrumental variable additive hazards models},
journal={Biometrics},
year={2015},
volume={71},
number={},
pages={122-130},
note={}
}
@ARTICLE{kjaersgaard,
author={Kjaersgaard, M. I. S. and Parner, E. T.},
title={Instrumental variable method for time-to-event data using a pseudo-observation approach},
journal={Biometrics},
year={2016},
volume={72},
number={},
pages={463-472},
note={}
}
@ARTICLE{mackenzie,
author={MacKenzie, T. A. and Tosteson, T. D. and Morden, N. E. and Stukel, T. A. and O'Malley, A. J.},
title={Using instrumental variables to estimate a Cox's proportional hazards regression subject to additive confounding},
journal={Health Services and Outcomes Research Methodology},
year={2014},
volume={14},
number={},
pages={54-68},
note={}
}
@ARTICLE{robins_tsiatis,
author={Robins, J. M. and Tsiatis, A. A.},
title={Correcting for non-compliance in randomized trials using rank preserving structural failure time models},
journal={Communications in Statistics - Theory and Methods},
year={1991},
volume={20},
number={},
pages={2609-2631},
note={}
}
@ARTICLE{loeys,
author={Loeys, T. and Goetghebeur, E.},
title={A causal proportional hazards estimator for the effect of treatment actually received in a randomized trial with all-or-nothing compliance},
journal={Biometrics},
year={2003},
volume={59},
number={},
pages={100-105},
note={}
}
@ARTICLE{yu_jrssb,
author={Yu, W. and Chen, K. and Sobel, M. E. and Ying, Z. L.},
title={Semiparametric transformation models for causal inference in time-to-event studies with all-or-nothing compliance},
journal={Journal of the Royal Statistical Society, Series B},
year={2015},
volume={77},
number={},
pages={397-415},
note={}
}
@ARTICLE{dabrowska,
author={Dabrowska, D. M. and Doksum, K. A.},
title={Partial likelihood in transformation models with censored data},
journal={Scandinavian Journal of Statistics},
year={1988},
volume={18},
number={},
pages={1-23},
note={}
}
@ARTICLE{price,
author={Price, A. L. and Patterson, N. J. and Plenge, R. M. and Weinblatt, M. E. and Shadick, N. A. and Reich, D.},
title={Principal components analysis corrects for stratification in genome-wide association studies},
journal={Nature Genetics},
year={2006},
volume={38},
number={},
pages={904-909},
note={}
}
@book{pft_book,
author = {Altalag, A. and Road, J. and Wilcox, P.},
title = {Pulmonary Function Tests in Clinical Practice},
publisher = {Springer-Verlag London},
year = 2009,
address = {London, United Kingdom},
note = {},
isbn = {}
}
@ARTICLE{terza_2008,
author={Terza, J. V. and Basu, A. Rathouz, P. J.},
title={Two-stage residual inclusion estimation: addressing endogeneity in health econometric modeling},
journal={Journal of Health Economics},
year={2008},
volume={27},
number={},
pages={531-543},
note={}
}
@ARTICLE{wooldridge_2015,
author={Wooldridge, J. M.},
title={Control function methods in applied econometrics},
journal={The Journal of Human Resources},
year={2015},
volume={50},
number={},
pages={420-445},
note={}
}
@ARTICLE{hanagiri,
author={Hanagiri, T. and Sugio, K. and Mizukami, M. and Ichiki, Y. and Sugaya, M. and Yasuda, M. and Takenoyama, M. and Yasumoto, K.},
title={Significance of smoking as a postoperative prognostic factor in patients with non-small cell lung cancer},
journal={Journal of Thoracic Oncology},
year={2008},
volume={3},
number={},
pages={1127-1132},
note={}
}
@ARTICLE{sardari,
author={Sardari Nia, P. and Weyler, J. and Colpaert, C. and Vermeulen, P. and Van Marck, E. and Van Schil, P.},
title={Prognostic value of smoking status in operated non-small cell lung cancer},
journal={Lung Cancer},
year={2005},
volume={47},
number={},
pages={351-359},
note={}
}
@ARTICLE{sawabata,
author={Sawabata, N. and Miyoshi, S. and Matsumura, A. and Ohta, M. and Maeda, H. and Sueki, H. and Hayakawa, M. and Okumura, M. and Sawa, Y.},
title={Prognosis of smokers following resection of pathological stage I non-small-cell lung carcinoma},
journal={General Thoracic and Cardiovascular Surgery},
year={2007},
volume={55},
number={},
pages={420-424},
note={}
}
@ARTICLE{gwas_catalog,
author={MacArthur, J. and Bowler, E. and Cerezo, M. and Gil, L. and Hall, P. and Hastings, E. and Junkins, H. and McMahon, A. and Milano, A. and Morales, J. and Pendlington, Z. M. and Welter, D. and Burdett, T. and Hindorff, L. and Flicek, P. and Cunningham, F. and Parkinson, H.},
title={{The new NHGRI-EBI Catalog of published genome-wide association studies (GWAS Catalog)}},
journal={Nucleic Acid Research},
year={2017},
volume={45},
number={},
pages={D896-D901},
note={}
}
@ARTICLE{huang_pnas,
author={Huang, Y. T. and Lin, X. and Liu, Y. and Chirieac, L. R. and McGovern, R. and Wain, J. and Heist, R. and Skaug, V. and Zienolddiny, S. and Haugen, A. and Su, L. and Fox, E. A. and Wong, K. K. and Christiani, D. C.},
title={Cigarette smoking increases copy number alterations in nonsmall-cell lung cancer},
journal={Proceedings of the National Academy of Sciences},
year={2011},
volume={108},
number={},
pages={16345-16350},
note={}
}
@ARTICLE{zhou_ccr,
author={Zhou, W. and Heist, R. S. and Liu, G. and Asomaning, K. and Miller, D. P. and Neuberg, D. S. and Wain, J. C. and Lynch, T. J. and Christiani, D. C.},
title={Second hand smoke exposure and survival in early-stage non-small-cell lung cancer patients},
journal={Clinical Cancer Research},
year={2006},
volume={12},
number={},
pages={7187-7193},
note={}
}
@ARTICLE{martinussen_2017,
author={Martinussen, T. and Vansteelandt, S. and Tchetgen Tchetgen, E. J. and Zucker, D. M.},
title={Instrumental variables estimation of exposure effects on a time-to-event endpoint using structural cumulative survival models},
journal={Biometrics},
year={2017},
volume={},
number={},
pages={},
note={DOI:10.1111/biom.12699}
}
@ARTICLE{stensrud_arx,
author={Stensrud, M. J.},
title={{Handling survival bias in proportional hazards models: A frailty approach}},
journal={arXiv},
year={2017},
volume={},
number={},
pages={},
note={arXiv:1701.06014}
}
@ARTICLE{stensrud_epi,
author={Stensrud, M. J. and Valberg, M. and Roysland, K. and Aalen, O. O.},
title={{Exploring selection bias by causal frailty models: The magnitude matters}},
journal={Epidemiology},
year={2017},
volume={28},
number={},
pages={379-386},
note={}
}
@ARTICLE{burgess_thompson,
author={Burgess, S. and Thompson, S. G.},
title={{Bias in causal estimates from Mendelian randomization studies with weak instrument}},
journal={Statistics in Medicine},
year={2011},
volume={30},
number={},
pages={1312-1323},
note={}
}
@ARTICLE{cao_genepi,
author={Cao, T. and Rajan, S. S. and Wei, P.},
title={Mendelian randomization analysis of a time-varying exposure for binary disease outcomes using functional data analysis methods},
journal={Genetic Epidemiology},
year={2016},
volume={40},
number={},
pages={744-755},
note={}
}
@ARTICLE{egger,
author={Bowden, J. and Davey Smith, G. and Burgess, S.},
title={Mendelian randomization with invalid instruments: effect estimation and bias detection through Egger regression},
journal={International Journal of Epidemiology},
year={2015},
volume={44},
number={2},
pages={512-525},
note={}
}
@ARTICLE{lin_ying,
author={Lin, D. Y. and Ying, Z.},
title={Semiparametric analysis of the additive risk model},
journal={Biometrika},
year={1994},
volume={81},
number={},
pages={61-71},
note={}
}
@ARTICLE{aalen_2015,
author={Aalen, O. O. and Cook, R. J. and Roysland, K.},
title={Does Cox analysis of a randomized survival study yield a causal treatment effect? },
journal={Lifetime Data Analysis},
year={2015},
volume={21},
number={4},
pages={579-593},
note={}
}
@ARTICLE{cho_huang,
author={Cho, S. H. and Huang, Y. T.},
title={Mediation analysis with causally ordered mediators using Cox proportional hazards model},
journal={Statistics in Medicine},
year={2019},
volume={},
number={},
pages={in press},
note={}
}
@ARTICLE{hernan_hr,
author={Hernan, M.},
title={The hazards of hazard ratios},
journal={Epidemiology},
year={2010},
volume={21},
number={},
pages={13-15},
note={}
}
@ARTICLE{,
author={},
title={},
journal={},
year={},
volume={},
number={},
pages={},
note={}
}
@ARTICLE{,
author={},
title={},
journal={},
year={},
volume={},
number={},
pages={},
note={}
}
PK×»‰ÐNÐNPK-
nTBWÒ}Šíöíö @MR_survival_sim_rev2.texUTñœePK-
nTBW×»‰ÐNÐN @@÷MR_survival.bibUTñœePK•ZF