Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
56,551 characters · 13 sections · 63 citation commands
Missing at Random or Not: A Semiparametric Testing Approach
\thispagestyle{empty}
Missing data are common in almost all data collection procedures where some units and/or items are not observed due to various reasons. How to appropriately deal with missing data has been a pervasive and crucial consideration for drawing valid scientific conclusions with statistical efficiency; see, among others, the monograph by Little2019.
In healthcare research, missing data are almost unavoidable, while their potential to undermine the validity of research results cannot be ignored. For example, the electronic health records (EHRs) data, as a result of not having been collected specifically for research purposes, are subject to considerable missing data. Missing data can be due to a lack of collection (e.g., patient was never asked about a condition) or a lack of documentation (e.g., patient was asked about a condition but the response was never recorded in the medical record). Lack of documentation is particularly common when it comes to a patient not having a symptom/comorbidity. Instead of recording a negative value for each potential symptom/comorbidity, all data fields are left blank (missing) and only the positive values are recorded. Thus it can be difficult to differentiate between the lack of a comorbidity, the lack of documentation of a comorbidity and the lack of data collection regarding the comorbidity. On the other hand, failure to account for the missing data can have a significant effect on the research conclusions sterne2009multiple,little2012prevention.
A central methodological concern in the missing data literature is on the mechanism governing the data missingness. There have been intensively developing investigations in this area, methodologically, theoretically, and practically. Incorporating different types of the data missingness mechanism is crucially relevant for developing appropriate statistical methods. In particular, the notions of missing completely at random, missing at random, and missing not at random are the commonly used classes for mechanisms characterizing the data missingness; see rubin1976 and Little2019. Missing completely at random refers to the case that events leading to any particular data being missing are independent of both observed and unobserved variables. Such a case is simple, but rarely happens in practice.
Missing at random supports a popular class of devices for investigating missing data problems. It conceptually refers to that the data missingness depends on observed variables alone but not the unobserved ones. A great deal of effort has been devoted in methodological development under this case; and many approaches have become routines for solving practical problems. Foremost, missing at random is convenient. For some cases in this scenario, applying existing methods with no adjustment by simply ignoring the missing data may lead to valid results; concerns are then on the efficiency and valid variance estimations for statistical inference; see, for example, White2010 and Bartlett2014. Missing at random is meritorious for providing unified statistical frameworks. That is, approaches of standard form can be developed and applied in broad class studies. For example, the likelihood approaches, Bayesian approaches, imputation methods are abundantly available; see the monographs by Molenberghs2007 and MFKTV_2014 for overview. For semiparametric methods using the generalized estimating functions LiangZeger_1986_Bioka, as another class of examples, if data are missing at random and the missing propensity function is appropriately modeled and estimated, the inverse probability weighting approach can ensure the validity of the resulting inference. For achieving the semiparametric estimation efficiency bound with missing data, one may design approaches by using the augmented inverse probability weighting approaches; see RRZ_1995_JASA, Tsiatis_2006, KimShao_2013 and references therein.
Missing not at random is most general; it refers to the case where the missingness is allowed to depend on the variables that are missing. Problems with missing not at random are more challenging. A main reason one may imagine is that there are plenty of individual possibilities associated with data missing not at random in various scenarios. Practically, these problems require dedicated effort and they are typically solved individually with relatively smaller class of situations compared with those under missing at random. As a consensus, when data are missing not at random, the propensity model plays a more important role, so that it is a fundamental part of the model; see the overview in monographs by Molenberghs2007 and MFKTV_2014.
For addressing challenging tasks, handling data missing not at random requires extra support from either stronger model assumptions and/or data structural information. For example, KottChang_2010 proposed some calibration approach in the context of survey sampling utilizing known totals of some variables; KimYu_2011_JASA investigated semiparametric regression with an exponential tilting modeling containing unknown parameters whose estimation requires extra information. It is even more challenging that these assumptions often cannot be validated without extra information and/or additional data collection procedure such as following-up studies.
Recently, a class of methods by using the instrumental variables are actively developing. Here we note that the instrumental variables are respecting to the data missing mechanism, whose objective differs in some ways from the conventional notion of instrumental variables as in econometrics and other areas where the data model is of the major concern. Among them, Wang2014 considered generalized methods of moments, and handled missing not at random with some estimating functions utilizing instrumental variables. ZhaoShao_2015_JASA proposed a novel approach for estimating generalized linear models under data missing not at random, without requiring the specification of a model for missing data mechanism. ShaoWang_2016 incorporated instrumental variables in inverse probability weighting, and relaxed requirement of extra information for the semiparametric approach of KimYu_2011_JASA. Riddles2016 developed a propensity score adjusted likelihood approach. Miao2016 investigated double robustness of the estimations with instrumental variable approaches. Morikawa2017 studied optimality of the estimation in a setting with instrumental variables. Zhao2018 constructed optimal pseudo-likelihood approach with data missing not at random. Wang2019 recently considered propensity selection and data modeling with some penalized information criteria. Tchetgen2018 studied identifiability in semiparametric estimation with instrumental variables; see also Miao2018.
Usually with missing not at random, correct specification of a model for the data missingness mechanism is required to ensure valid and/or efficient inference. However, such a specification is generally difficult as it closely depends on the part of the data that are not observed. In addition, methods to address missing not at random can appreciably increase the complexity and the uncertainty of study results; so that its validity is clearly more desirable. Nevertheless, relative to the recent surge of the development in studying missing not at random, few methodology is available for conducting statistical testing against the propensity model specifications. Mohan2014 considered testing for data missing mechanism in the context of directed acyclic graph and causal relationships between variables; they established conditions when such testing is possible. Recently, Breunig2019 proposed a method with instrumental variables for testing the missing at random assumption with some squared integrated distance built upon some conditional moment identities.
We consider in this study a new approach for objectively deciding the mechanism of the data missingness from the two classes: missing at random or missing not at random. Such an approach has practical interest to an applied scientist who can be reluctant, without concrete supporting evidence, to undertake analysis methods for missing not at random. Inspired by recent development in handling missing not at random with instrumental variables, we develop a semiparametric hypothesis testing framework with this class of methods. To stay focused in our presentation, we consider the generalized linear models GLM_1989_MN, while recognizing the principle of the development broadly applies. For the mechanism governing the data missingness, we are motivated by the setting of ZhaoShao_2015_JASA where the specification of a model for the missing data mechanism is not required, entitling broad validity of the resulting estimator. Then inspired by the rationale of the famous Hausman's test hausman1978specification, we propose to identify two parameter estimators that both are valid if missing at random holds, and only one is valid otherwise. Then the signal for our test relies on a discrepancy measure taking significantly larger values when missing not at random is more appropriate. There are many candidates for the first estimator, and in our study we take the inverse probability weighted estimator as a concrete example in our development; and we observe that our method broadly applies. For the second estimator, we apply the semiparametric approach of ZhaoShao_2015_JASA that is valid for both classes of data missingness. Our theory confirms that the testing procedure is valid and powerful, and our simulation studies show that the test works satisfactorily when assuming data missing at random is not appropriate.
Our investigation makes a few contributions. Practically, we provides a framework that missing at random assumption can be tested against a broad setting. Then to what end such an assumption is reasonable can actually be evaluated. Methodologically, our approach attempts to adequately exploits the current setting for data missing not at random, and we demonstrate that testing the missing at random assumption is possible. As for technical development, existing Hausman's tests are developed with full parametric models, while our development extends their applicability to broad semiparametric settings. We analytically calculate the influence function of the semiparametric estimator of ZhaoShao_2015_JASA, which provide a tool for broad statistical inferences. Our development accommodating a nonparametric component as an infinite dimensional nuisance parameter is a new feature of its own interest in the context of model specification tests in the context of missing data problems.
The rest of this article is organized as follows. In Section 2, we propose a general testing approach and provide a specific test statistic as an example. We showed the limiting distribution of the test statistic is a $\chi^2$ distribution. In Section 3, we apply the proposed test to data from the Keep-it-off study to investigate the missing mechanism of the self-reported body weight of participants in a weight loss program. In Section 4, we present a comprehensive simulation study to show the performance of the proposed test in terms of type I error and power. We provide a discussion in Section 5, and we include the regularity conditions and part of the proofs in the Appendix.
We concretely consider a generic setting of regression analysis in this study. For a data set with a response variable $Y$ and covariates $X = (U^{T}, Z^{T})^{T}$ (the difference between $U$ and $Z$ to be explained later), we consider a setting with generalized linear model. Specifically, the conditional distribution of $Y$ given covariates $X$ belongs to the exponential dispersion family, and the conditional mean of $Y$, $\mu=E (Y|X)$, is related to covariates through a link function $h(\cdot)$, i.e.
where $\eta$ is the natural parameter, $\lambda$ is the dispersion parameter, $\beta = (\alpha, \beta_u^{T}, \beta_z^{T})^{T}$ are regression coefficients, and $\mu(\eta)=b'(\eta)$ by the property of the exponential family. We denote the dimensions of $U$, $Z$ and $X$ by $m_u$, $m_z$ and $m$ respectively.
As in econometrics literature Heckman1979,Heckman1997, instrumental variables originally refer to those satisfying two conditions: 1) they are conditional independent of the outcome variables that can be missing, given other variables; and 2) they are affecting the missing mechanism in a model with all variables. In recent missing data literature, conditional independence is assumed between the instrumental variables and the missingness given all other variables including the missing ones. Additionally, some identifiability conditions are needed for the parametric data model on the connections between the instrumental variables and the response variables; see, for example, Fang2016 and Wang2019.
We consider that the response variable $Y$ can be missing, and let $R_i=1$ if $Y_i$ is observed, and $0$ otherwise. We consider the following two assumptions as the possible underlying data missingness mechanisms. \\ {Assumption A1:} \quad $P(R|Y, U, Z)=P(R|U)$;\\ {Assumption A2:} \quad $P(R|Y, U, Z)=P(R|Y, U)$.
Assumption A1 refers to the case of missing at random, where the missingness is conditional independent of $Y$ that can be missing. Assumption A2 corresponds to the case of missing not at random. Here the component $Z$ is referred to as the instrumental variables. Clearly, A1 is stronger, and itis a special case of A2. Thus, we note that a valid estimator under A2 remains valid under A1. Here common to both assumptions is that, conditioning on $Y$ and $U$, the variable $Z$ is independent of the missing data indicator $R$.
Estimating model parameters with missing at random has been intensively studied, and commonly applied approaches include the inverse probability weighting, augmented inverse probability weighting, multiple imputations and others; see Tsiatis_2006, Molenberghs2007, and Little2019. When data are missing not at random, inference on $\beta$ commonly requires specification of the missing propensity function. Nevertheless, with the instrumental variable $Z$, valid estimation of the parameter $\beta$ can be made semiparametrically without a need to specify a model for the missingness propensity function ZhaoShao_2015_JASA, subject to some identifiability conditions.
We now define some notations. For any function $t(\theta,X,Y,R)$, we use $\nabla_\theta t(\theta,X,Y,R) = \partial t(\theta,X,Y,R)/\partial \theta$, and $\nabla_{\theta\theta} t(\theta,X,Y,R) = \partial^2 t(\theta,X,Y,R)/\partial \theta^2$. We understand $t_i$ as $t(\theta,X,Y,R)$ with $(X,Y,R)$ replaced by $(x_i,y_i,r_i)$, and treat $t_i$ as a random variable $(i=1,\dots,n)$.
In practice, it is important to understand what the true underlying missing data mechanism is, which can be formulated as testing Assumption A1 against Assumption A2. Since A1 is a special case of A2, a rationale to develop a statistical test is then identifying two estimators for $\beta$, $\widehat \beta$ and $\widetilde \beta$, such that $\widehat\beta$ is valid only under assumption A1, and $\widetilde\beta$ is valid under assumption A2. Clearly, $\widehat\beta$ is expected to be biased when A1 is violated; and $\widetilde\beta$ is valid under both A1 and A2. Subsequently, the discrepancy between $\widehat \beta$ and $\widetilde \beta$ becomes the signal for detecting violation of Assumption A1. When such discrepancy is large, we have evidence that missing at random is suspicious. Such a rationale is the spirit of the classical Hausman’s test hausman1978specification, though the test was developed for parametric estimators and did not formally address estimators that require nonparametric estimation of nuisance parameters.
In our framework, a main challenge is that we have a broad assumption of the form A2. Here we intend to impose minimal restriction on the form of the propensity function -- allowing it to be nonparametric -- so that our approach is broadly applicable. It is worth to note that that conventional Hausman's test does not apply due to this scope of our work incorporating a setting requires handling a nonparametric distribution. Indeed, as shown in our later development, profiling out this nonparametric component as a functional nuisance parameter involves major technical challenges for constructing the testing statistic. Therefore, as an interest of its own, our test statistic extends the scope of the classical Hausman's test with nonparametric nuisance parameters that can be infinite dimensional.
Concretely, we first show the following theorem describing the asymptotic property of the discrepancy measure $(\widehat\beta-\widetilde\beta)$ in a general setting accommodating semiparametric estimations. We denote by $\Lambda_1$ and $\Lambda_2$ nuisance parameters respectively under Assumption A1 and A2 of the specific approaches for estimating the model parameter $\beta$. Since model estimation approaches are broad including those parametric and semiparametric ones, here $\Lambda_1$ and $\Lambda_2$ are allowed to be general including functional-valued quantities such as nonparametric propensity functions and distribution functions in some semiparametric models. They essentially depend on the model settings under which the estimation methods are developed; see our examples in Section (ref) and (ref).
In Theorem (ref), $\psi(\cdot)$ and $\phi(\cdot)$ are known as the influence functions. Theorem (ref) ensures that to construct a test statistic, one needs to find a consistent estimator for the variance matrix $W$.
In existing Hausman's tests, if $\widehat{\beta}$ is efficient using the likelihood approach under Assumption A1, the variance matrix $W$ has the property that $W = V_2-V_1$, where $V_2$ is the variance of $\widetilde{\beta}$ and $V_1$ is the variance of $\widehat{\beta}$. In such a case, the application of the Hausman's test is straightforward.
More generally, however, there are two major difficulties in our study. First, choices of the estimators are broad, and semiparametric approach are popular in existing methods. Thus, calculating the variances of $\widehat \beta$ and $\widetilde \beta$ is generally involved due to nonparametric nuisance parameters. Second, our framework accommodates cases when $\widehat \beta$ is not necessarily efficient. Thus $\widehat \beta$ and $\widetilde \beta-\widehat \beta$ could be correlated, and it no longer holds with a simple form that $W = V_2-V_1$. To tackle this challenge, we propose to apply the sample covariance matrix of the difference between the influence functions, i.e., $\psi(y_i,x_i; \beta_0, \Lambda_1)-\phi(y_i,x_i; \beta_0, \Lambda_2)$, $(i=1,\dots, n)$ as the estimator for $W$. With a consistent estimator $\widehat W$, we can show that the test statistic $T = n(\widehat\beta-\widetilde\beta){\widehat W}^{-1}(\widehat\beta-\widetilde\beta)$ converges in distribution to $\chi^2_{m+1}$.
There are broad choices for estimating the model parameter when data are missing at random. In this case, even the complete-case analysis ignoring all missing data may provide a valid inference White2010. To make our method more generalizable, we consider here an exemplary method for missing at random data, which is the inverse probability weighting (IPW) method with nonparametric estimation of the propensity function.
When the true propensity $\pi_0(u) = P(R = 1|U = u)$ is known, the inverse probability weighted likelihood score equation can be constructed as
where $S(Y,X;\beta)= {\nabla_{\beta}\log p(Y|X;\beta)}$. However, in practice, the propensity function is unknown. In this case, parametric models can be specified for the missing propensity, and a two-stage method can be used to obtain a plug-in estimating equation. To avoid misspecification of the propensity function, we consider estimating $\pi_i$ through nonparametric regression. Let $J(x)$ be a multivariate kernel function of variable $u$ satisfying $\int J(u) du = 1$. The missing propensity $\pi_0(u)$ can be estimated by a kernel regression estimator of form
where $b$ is the bandwidth that determines smoothness of $\widehat{\pi}_b$. By plugging in the estimator $\widehat {\pi}_b(u)$ back to estimating equation ((ref)), we obtain the following estimating equation
With a proper choice of bandwidth $b$, it can be shown that the kernel regression estimator $\widehat {\pi}_b(u)$ is a consistent estimator of $\pi_0(u)$. Although the convergence rate of $\widehat {\pi}_b(u)$ is slower than ${n}^{1/2}$ hardle1993comparing, as suggested by Newey_1994 and NeweyMcFadden1994 the estimator obtained from the plug-in estimating equation is still ${n}^{1/2}$-consistent and asymptotically normal. Specifically, we have the following theorem.
The proof of Theorem (ref) is outlined in the Supplementary Material. This is an example of parameter estimation with nonparametric propensity-weighted score functions. Here, it is the key to develop the influence function $\psi_i$ with the propensity function $\pi_0(u)$ estimated nonparametrically. Our methodological framework here also more broadly applies.
The missing not at random case is a more challenging scenario from multiple aspects; see, among others, KimShao_2013. Our attempt here has a consideration to be most accommodative with no extra assumption beyond A2. We have such a candidate, as elaborated below. As pointed out by ZhaoShao_2015_JASA, a pseudolikelihood approach can be built with the instrumental variables $Z$ that satisfy Assumption A2. A key observation is the following equivalence between the conditional distributions
where the first equation is from A2, and the second equation is simply its definition. Hence, from ((ref)) evaluating the likelihood is valid with complete only data such that $R=1$.
To avoid a fully parametric specification of the model, the joint distribution $P(U,Z)$ becomes the nonparametric nuisance parameter -- a major difference from the development in Section (ref) with missing at random. In this context, we consider the nonparametric product kernel estimator
where $K$ is a one-dimensional kernel function, $h$ is the bandwidth and $K_m$ denotes the $m$-dimensional kernel function. A plug-in semiparametric log-pseudolikelihood function is then obtained
where $\widehat{F}$ is the corresponding empirical cumulative distribution function of $\widehat f$. Then the maximum pseudolikelihood estimator is $\widetilde{\beta}={\textrm{argmax}_{\beta}} \left\{ \ell(\beta,\widehat{F}) \right\}$. We also have that the first order derivative function \[ h(\beta,\widehat{F}) = \frac{1}{n} \sum_{i=1}^{n}{\nabla_\beta H_i(\beta,\widehat{F})} =\frac{1}{n} \sum_{i=1}^{n}r_i\left\{ \frac{\nabla_{\beta}p(y_i|u_i,z_i;\beta)}{p(y_i|u_i,z_i;\beta)} -\frac{\int \nabla_{\beta}p(y_i|u_i,z;\beta)\widehat{f}(u_i,z)dz}{\int p(y_i|u_i,z;\beta)\widehat{f}(u_i,z)dz}\right\}, \] and let $\mu(F) = E\{\nabla_\beta H(\beta_0,F)\}$ for a given distributional function $F$ of $U$ and $Z$.
To construct a discrepancy measure between two estimators $\widehat \beta$ and $\widetilde \beta$, both of their influence functions are required. The influence function of $\widetilde \beta$ is more complicated, and its explicit form is not known in the literature. To solve this difficulty, we establish the following lemma, providing the explicit form of the approximation error between the nonparametric estimation and the truth, which is a key component in the influence function of $\widetilde \beta$.
The proof of Lemma (ref) is outlined in Appendix B where appropriately handling the nonparametric joint distribution $P(U,Z)$ is the key. As an intermediate result for constructing our test, we present the following theorem, refining the results of ZhaoShao_2015_JASA by providing the asymptotic distribution of $\widetilde{\beta}$, as well as an explicit form of its influence function.
A proof is given in Appendix B. The form of $\phi_i$ in Theorem (ref) has clear insights. The contribution in $\phi_i$ without $\delta_i$ reflects the accuracy of the oracle estimator as if the true distribution $F$ is known. The part from $\delta_i$ reflects the approximation error due to the nonparametric estimation. Although the nonparametric distribution estimation has a slower convergence rate than $n^{1/2}$, the parameter estimation still achieves the $n^{1/2}$ rate, echoing the phenomenon known in the literature on semiparametric estimations; see, for example Newey_1994 and NeweyMcFadden1994.
With the two estimators obtained in Sections 2.3 and 2.4, we construct by evaluating the discrepancy between $\widehat \beta$ and $\widetilde \beta$ such that the test statistic is expected to take a large value when Assumption A1 is violated. Concretely, let $T$ be a test statistic defined as
where the appropriate weighting matrix $W$ is constructed from the estimated influence functions $\phi_i$ and $\psi_i$ as in Theorems (ref) and (ref):
The estimators of $\phi_i$ and $\psi_i$ are respectively defined as
and
with
The key property of the test statistic $T$ in equation ((ref)) is described in the following corollary.
A proof is provided in the Supplementary Material. Corollary (ref) reveals a central $\chi^2$ limiting distribution of $T$ under the null hypothesis, making it convenient for practical implementation of the proposed test.
We then evaluate the asymptotic power of the test statistic $T$ in equation ((ref)) under a broad class of local alternatives. Suppose the missing propensity relates to $Y$ through a parameter $\gamma$, which is with finite dimension $q\ge 1$. We denote the propensity by $P (R=1|y, u) = \pi(y,u;\gamma)$. Suppose when $\gamma = 0$, the missingness is at random and $\pi(y,u;\gamma)$ reduces to $\pi_0(u)$. Assume that $\pi(Y,U;\gamma)$ is second order differentiable, and there exist $\epsilon>0$ such that $ \pi(y,u;\gamma) \in (\epsilon,1-\epsilon)$ for all $\gamma$. Consider a fixed $\gamma_0$ and a sequence of alternatives:
With LeCam's third lemma Van_1998, we can establish the following results.
A proof of the theorem is given in the Supplementary Material. The parameter $\gamma$ represents the dependence of the missingness propensity function on $Y$. Theorem (ref) implies that the test statistic remains bounded under the alternative $H_\alpha:\gamma = n^{-1/2}\gamma_0$, and its power is determined by the non-centrality parameter $\gamma_0^{T}\eta^{T}W^{-1}\eta\gamma_0$. Furthermore, when the signal of missing not at random is not vanishing to zero at a rate faster than $n^{1/2}$, the power of the test with $T$ will go to one.
We analyzed data from the Keep It Off randomized controlled trial which recruited 191 participants and collected their body weight over a 1-year period shaw2017design. Participants were randomized to receive one of three interventions to help maintain their weight after a recent weight loss. The primary outcome was the participant's weight at the in-person weigh-in at month six. As part of the follow-up, the participants reported their daily body weights through use of a wireless scale that transmitted data to the Way to Health portal, an online clinical study platform WTH2009. For more detailed information about the design of the trial, we refer to shaw2017design. In this example, we investigated the nature of the missing data for the daily at-home weigh-in data for the combined cohort, with the help of a regression analysis between daily weight and baseline covariate variables collected in the trial.
We chose the self-reporting at-home weights at the day before the month-six visit as the outcome ($Y$). The missing proportion of the at-home weights is around 47%. Since the primary outcome, the participant's weight at the in-person weigh-in at month six, was measured for all participants, we considered it as a gold standard ($Y^\star$) since it is reasonable to assume no large variation of weight within two days. The availability of the gold standard $Y^*$ in this scenarios allows us to validate our methods.
We observed variables including age, gender, race, education-level, baseline weight, baseline boday mass index (BMI), treatment arms, and self-reported scores of physical activity and eating habits based on a questionnaire provided by the Way to Health portal. In the existence of the gold standard in person weight at month six, we first fitted a linear regression using all the aforementioned variables as covariates. Among all these variables, only baseline weight and baseline BMI are significantly associated with weight at month six, and this finding is also consistent with many existing findings ben2003predictors. Hence, to reduce the complexity of the proposed test, we chose to use the linear regression between daily weight and the baseline weight and BMI.
Among baseline weight ($U$) and BMI ($Z$), baseline BMI was hypothesized to be an instrumental variable, which should be (1) associated with the daily weight, and (2) after conditioning on the daily weight and the baseline weight, the missingness of the daily weight at a specific day should not be related to the baseline body mass index. To investigate hypothesis (1), we fitted a linear regression between the gold standard $Y^*$ and $X = (U,Z)$ and the results were shown in Table (ref). We observed that the body mass is significantly associated with the outcome $Y^*$ (p-value $<$ 0.01).
To validate assumption (2), we fitted the following logistic regression
where $R = 1$ if $Y$ is observed. Specifically, Table (ref) below showed the model fitting results of the above logistic regression, where the coefficient of the potential instrumental variable, baseline body mass index ($Z$), is statistically non-significant (p=0.26). The observation is reasonable that the current weight and baseline weight will capture information regarding current weight loss, likely a dominant determinant in whether a participant decides to weigh-in on a given day.
Applying our method using the baseline body mass index as an instrumental variable led to a test statistic of 9.12, compared to a chi-squared distribution with $3$ degrees of freedom, corresponding to a p-value of 0.02. This test suggested that the data missingness is not at random. This matched the observation from the data that people having smaller body weight are more likely to report their weight for a given body mass index and baseline weight. This is consistent with the conventional wisdom that participants in weight loss trials that are gaining weight in the trial are disappointed at the lack of success and more likely to skip reporting.
This conclusion is also validated by the result in Table (ref). Specifically, we found that conditioning on baseline weight and body mass index, weight at month 6 ($Y^\star$), as an approximation to $Y$, remains statistically significantly associated with the missing indicator (p-value $=$ 0.02), which also suggested that the data are not missing at random.
Finally, Figure 1 provided a graphical evidence of the missing not at random assumption. In the left panel, the missing propensity estimated from model ((ref)) was plotted against the missing propensity estimated from the following model assuming missing at random.
Substantial differences of estimated propensities under the two assumptions were observed. The largest difference was as high as 0.21, and 24% of the differences were observed to be greater than 0.1. In the right panel, we plotted the log odds ratio of the estimated propensity from model ((ref)) against the estimated propensity from model ((ref)) versus the gold standard weight. From the figure, we observed a clear decreasing trend of the log odds with the increasing of the true body weight, which implied the difference of the two estimated propensities were related to the true body weight. For participant with less weight ($< 200$ lb), model ((ref)) with missing at random tended to overestimate the propensity, and for larger weight, model ((ref)) tended to underestimate the propensity. These results were consistent with the conclusion of missing not at random obtained by both our proposed test and the logistic regression in Table 2.
We conducted simulation studies to examine the finite sample behavior of our proposed test. We considered univariate $Y$, $U$ and $Z$ from multivariate normal distribution where $[Y|U,Z] \sim N(1+U+b_zZ,1)$, $[U|Z] \sim N(1-Z,1)$ and $Z\sim N(0,1)$. We generated the missingness indicator variable $R$ for the response variable $Y$ from a Bernoulli distribution with $pr(R=1|Y,U,Z)\sim \Phi(c_0+c_1f(Y)+c_2U)$, where $\Phi(\cdot)$ is the cumulative density for standard normal distribution. By specifying different functional forms for $f(\cdot)$, we allowed different types of dependence of missingness on the value of $Y$. The parameter $b_z$ can be considered as a quantification of how strong the instrumental variable is related to $Y$. In the simulation, we let $b_z = 1$ or $0.5$ to investigate the impact of association between $Z$ and $Y$ on the power of the proposed test. Furthermore, the magnitude of $c_1$ specifies the strength of dependence of missingness on $Y$. The value $c_1=0$ corresponds to the null hypothesis setting where the missingness is at random. We varied the value of $c_1$ to evaluate the size and power for the proposed test. We also chose different values for $c_2$ in order to evaluate the influence of the strength of association between missingness and the variable $U$ on the size and power of the test. An appropriate value of $c_0$ was chosen for each $(c_1, c_2)$ to get an overall missing percentage of $20\%$. Details of the parameter settings are proposed in Table S1 of the Supplementary Material.
Table (ref) showed the estimated levels of Type I error and power from $1000$ replications of the test with sample size of 1000. When the null hypothesis of missing at random is true, i.e. $c_1=0$, most of the rejection rates of the test were within $95\%$ confidence intervals for the nominal levels of $5\%$, i.e. $3.6$--$6.4\%$. A plot of the quantiles of the test statistics against those of the chi-squared distribution indicated that the chi-squared distribution works well at levels other than $5\%$; see Figure S1 in the Supplementary Material. In Table 4, we observed that when sample size increased to $2000$, the test still controlled Type I error at nominal level and obtained higher statistical power.
The power of the test of missing at random increased with $c_1$, which quantifies the dependence between missingness and response variable. The functional form $f(\cdot)$ has a sizable impact on the power of the proposed test. On the other hand, the magnitude of the parameter $c_2$ has relatively small impact on the power. At the nominal level of $5\%$, the proposed test had about $90\%$ power when $c_1$ was $0.3$ for $f(y)=y$, when $c_1$ was $0.4$ for $f(y)=0.4y^2$, and when $c_1$ was $0.4$ for $f(y)=2.5I(y>1)$. The magnitude of the parameter $b_z$ has a strong impact on the power of the proposed test. When $b_z$ reduced from $1$ to $0.5$, the power also dropped nearly a half for all scenarios.
In summary, the proposed test controlled the Type I error well and was reasonably powerful in detecting missing not at random in the settings we considered.
In this paper we proposed a new testing procedure to investigate the missing at random assumption for an outcome $Y$ in generalized linear models. We developed a general test statistic in Theorem 1 based on a discrepancy measure of two estimators under two different assumptions, respectively, i.e., missing at random and missing not at random. We provided a realization of the proposed test statistic by choosing the IPW estimator under the missing at random assumption and the estimator proposed by ZhaoShao_2015_JASA under the missing not at random assumption constructed {with} the existence of an instrumental variable. Using a newly developed method by ichimura2015influence, we derived the influence functions of the two estimators in a semiparametric setting to avoid parametric specification of the missing propensity and the joint distribution of the covariate variables. The realization of the proposed test was validated and evaluated by a simulation study and we found the proposed test was able to control the type I error and provide reasonable statistical power. {Using data from a weight loss study where 47% of the at-home body weights of 191 participants were missing}, we investigated the nature of the missing data mechanism using BMI as an instrumental variable. We found strong evidence of missing not at random. Such a finding is consistent with our analysis using validation data where the body weights of all participants were measured a day after the at-home body weights. Our cases study illustrated the practical utility of the proposed testing procedure. The R code has been properly documented and is available online.
The general testing procedure in Theorem 1 covers a broad class of {tests for investigating the missing data mechanism} and can be extended beyond the scope of the generalized linear models. For example, we can consider the semiparametric density ratio model by luo2011proportional which extends the generalized linear model by leaving the reference distribution unspecified. When the outcome is longitudinal, models proposed in luo2014moment and chen2015regression can also be considered. While all being flexible, a challenge is to find two estimators and derive their corresponding influence functions, where one is only valid under missing at random and the other is valid under both missing at random and missing not at random assumptions.
In this paper, the dispersion parameter $\lambda$ is considered as known in our paper. If $\lambda$ is unknown, we define $\theta = (\beta^{T},\lambda)^{T}$, and all the conclusions in this paper still hold with $\beta$ replaced by $\theta$ with corresponding adjustments in the regularity conditions. Alternatively, we may also adapt a two-stage estimation procedure via pseudolikelihood gong1981pseudo,liang1996asymptotic,chen2010asymptotic, where $\lambda$ can be replaced by a consistent estimator. In addition, from the perspective of estimation rather than hypothesis testing, the quantity $(\widetilde \beta - \widehat \beta)/{\widehat \beta}$ can be used to measure the relative bias attributable to the assumption of missing at random when the missingness is truly not at random. This measure can be a helpful complement to a p-value from the proposed testing procedure. Finally, the proposed procedure can be extended to the scenarios with missing covariate variables. Some of the extensions are currently under investigation and will be reported in the future.
\setstretch{1.24}