Semiparametric Analysis of the Proportional Likelihood Ratio Model and Omnibus Estimation Procedure

TL;DR

This paper conducts a semi-parametric analysis of the proportional likelihood ratio model, deriving efficient score functions and proposing estimators applicable to various missing-data scenarios, supported by simulation results.

Contribution

It introduces a family of Z-estimators, including an efficient one, for the proportional likelihood ratio model, with applicability to missing data mechanisms.

Findings

Proposed estimators perform well in finite samples.

Explicit efficient score functions derived for the model.

Estimators applicable to censored and non-random sampling data.

Abstract

We provide a semi-parametric analysis for the proportional likelihood ratio model, proposed by Luo & Tsai (2012). We study the tangent spaces for both the parameter of interest and the nuisance parameter, and obtain an explicit expression for the efficient score function. We propose a family of Z-estimators based on the score functions, including an approximated efficient estimator. Using inverse probability weighting, the proposed estimators can also be applied to different missing-data mechanisms, such as right censored data and non-random sampling. A simulation study that illustrates the finite-sample performance of the estimators is presented.

Tables

Table 1. Table 1: Simulation results with β 0 = ( − 1 , − 1 ) T subscript 𝛽 0 superscript 1 1 𝑇 \beta_{0}=(-1,-1)^{T} : empirical mean (empirical standard deviation)

	Setting 3
		100	200	400	800
${\beta }_1$	Lou & Tsai	-1.09 (0.94)	-0.99 (0.65)	-1.00 (0.43)	-1.00 (0.29)
	Chan	-1.10 (0.97)	-1.00 (0.67)	-1.01 (0.44)	-1.01 (0.30)
	Proposed	-1.07 (0.95)	-0.98 (0.66)	-1.00 (0.44)	-0.99 (0.30)
${\beta }_2$	Lou & Tsai	-1.06 (0.90)	-1.05 (0.61)	-1.02 (0.41)	-1.02 (0.27)
	Chan	-1.08 (0.92)	-1.06 (0.63)	-1.02 (0.42)	-1.02 (0.27)
	Proposed	-1.06 (0.90)	-1.05 (0.62)	-1.02 (0.42)	-1.02 (0.27)
Distance	Lou & Tsai	0.08 (0.04)	0.06 (0.03)	0.06 (0.02)	0.03 (0.01)
	Chan	0.08 (0.04)	0.06 (0.03)	0.04 (0.02)	0.03 (0.01)
	Proposed	0.08 (0.04)	0.06 (0.03)	0.04 (0.02)	0.03 (0.01)
	Setting 4
${\beta }_1$	Lou & Tsai	-1.06 (1.42)	-1.01 (0.93)	-1.00 (0.62)	-1.01 (0.43)
	Chan	-1.08 (1.47)	-1.03 (0.96)	-1.00 (0.64)	-1.02 (0.44)
	Proposed	-1.32 (4.32)	-1.21 (1.20)	-1.09 (0.78)	-1.05 (0.53)
${\beta }_2$	Lou & Tsai	-1.61 (1.30)	-1.56 (0.86)	-1.49 (0.58)	-1.47 (0.38)
	Chan	-1.65 (1.35)	-1.58 (0.89)	-1.50 (0.59)	-1.48 (0.39)
	Proposed	-1.31 (1.33)	-1.16 (0.93)	-1.01 (0.61)	-0.94 (0.41)
Distance	Lou & Tsai	0.13 (0.06)	0.12 (0.05)	0.11 (0.04)	0.11 (0.03)
	Chan	0.16 (0.06)	0.12 (0.05)	0.11 (0.04)	0.11 (0.03)
	Proposed	0.12 (0.01)	0.09 (0.07)	0.07 (0.05)	0.05 (0.03)

Table 2. Table 2: Simulation results with β 0 = ( − 1 , − 1 ) T subscript 𝛽 0 superscript 1 1 𝑇 \beta_{0}=(-1,-1)^{T} : empirical mean (empirical standard deviation)

	Setting 1
		100	200	400	800
${\beta }_1$	Lou & Tsai	-1.06 (0.77)	-1.00 (0.55)	-1 (0.37)	-1.01 (0.25)
	Chan	-1.07 (0.78)	-1.01 (0.56)	-1.00 (0.38)	-1.01 (0.25)
	Proposed	-1.06 (0.77)	-1.00 (0.55)	-1.00 (0.37)	-1.01 (0.25)
${\beta }_2$	Lou & Tsai	-1.05 (0.71)	-1.03 (0.51)	-1.01 (0.33)	-1.01 (0.22)
	Chan	-1.07 (0.72)	-1.03 (0.53)	-1.02 (0.34)	-1.01 (0.23)
	Proposed	-1.05 (0.71)	-1.03 (0.51)	-1.01 (0.33)	-1.01 (0.22)
Distance	Lou & Tsai	0.07 (0.04)	0.05 (0.03)	0.03 (0.02)	0.02 (0.01)
	Chan	0.07 (0.04)	0.05 (0.03)	0.03 (0.02)	0.02 (0.01)
	Proposed	0.07 (0.04)	0.05 (0.03)	0.03 (0.02)	0.02 (0.01)
	Setting 2
		100	200	400	800
${\beta }_1$	Lou & Tsai	-1.06 (0.24)	-1.03 (0.16)	-1.02 (0.11)	-1.01 (0.08)
	Chan	-1.05 (0.26)	-1.02 (0.17)	-1.01 (0.12)	-1.01 (0.08)
	Proposed	-1.06 (0.24)	-1.03 (0.16)	-1.02 (0.11)	-1.01 (0.08)
${\beta }_2$	Lou & Tsai	-1.07 (0.23)	-1.04 (0.16)	-1.02 (0.10)	-1.01 (0.08)
	Chan	-1.05 (0.25)	-1.03 (0.16)	-1.01 (0.11)	-1.01 (0.08)
	Proposed	-1.07 (0.23)	-1.04 (0.16)	-1.02 (0.10)	-1.01 (0.08)
Distance	Lou & Tsai	0.39 (0.20)	0.28 (0.14)	0.20 (0.10)	0.13 (0.07)
	Chan	0.40 (0.21)	0.29 (0.14)	0.20 (0.10)	0.14 (0.07)
	Proposed	0.39 (0.20)	0.28 (0.14)	0.20 (0.10)	0.13 (0.07)

Equations

$$p_{Y,X}(y,x)=p_{Y\mid X}(y\mid x)p_X(x)=\frac{\exp ({\beta }^{T}xy)g(y)}{\int \exp ({\beta }^{T}xy)dG(y)}\eta (x),$$

$$p_{Y\mid X}(y\mid x)=\frac{\rho (y,y_0,x,x_0,\beta )g(y)}{\int \rho (y,y_0,x,x_0,\beta )dG(y)},$$

$$\rho (y,y_0,x,x_0,\beta )=\frac{p_{Y\mid X}(y\mid x)}{p_{Y\mid X}(y\mid x_0)}{\left\{\frac{p_{Y\mid X}(y_0\mid x)}{p_{Y\mid X}(y_0\mid x_0)}\right\}}^{-1},$$

$$\mathcal{P}=$$

$$p_Y(y)=g(y)\int \frac{\exp ({\beta }^{T}xy)\eta (x)}{\int \exp ({\beta }^{T}xy)dG(y)}dx,$$

$$\text{(i)}{\Lambda }_1=\left\{h(Y)-E\{h(Y)\mid X\}:h(Y)\text{ is a q-dimensional vector-valued function}\right\}$$

$$\text{(ii)}{\Lambda }_2=\left\{\alpha (X):\alpha (X)\text{ is a q-dimensional vector-valued function such that }E\{\alpha (X)\}=0\right\}$$

$$E\left([h(Y)-E\{h(Y)\mid X\}]\alpha (X)\right)=$$

$$=$$

$$\mathcal{K}\{a(Y,X)\}=E\{a(Y,X)\mid Y\}-E[E\{a(Y,X)\mid X\}\mid Y].$$

$$(i)\Pi \{a(Y,X)\mid {\Lambda }_1\}={\left(I-\mathcal{B}\right)}^{-1}\mathcal{K}\{a(Y,X)\}-E\left[{\left(I-\mathcal{B}\right)}^{-1}\mathcal{K}\{a(Y,X)\}\mid X\right],$$

$$(ii)\Pi \{a(Y,X)\mid {\Lambda }_2\}=E\{a(Y,X)\mid X\}.$$

$$\Pi \{a(Y,X)\mid {\Lambda }^{\perp }\}=$$

$$-\sum _{j=0}^{\infty }\left({\mathcal{B}}^j\mathcal{K}\{a(Y,X)\}-E\left[{\mathcal{B}}^j\mathcal{K}\{a(Y,X)\}\mid X\right]\right).$$

$$XY-E(XY\mid X)-\sum _{j=0}^{\infty }\left[{\mathcal{B}}^j\mathcal{K}(XY)-E\{{\mathcal{B}}^j\mathcal{K}(XY)\mid X\}\right].$$

$$XY-E(XY\mid X)-E(XY\mid Y)+E\{E(XY\mid Y)\mid X\}.$$

$$\Pi \{\phi (Y,X)\mid \Lambda \}\equiv 0^{q\times 1}.$$

$$m_{0,a}(y,x)$$

$$=a(x,y)-E\{a(X,Y)\mid X=x\}-E\{a(X,Y)\mid Y=y\}$$

$$+E[E\{a(X,Y)\mid Y\}\mid X=x].$$

$$\begin{array}{rl}{S}_a(x,\beta ,G)& =\int a(s,x)\exp ({\beta }^{T}xs)dG(s),\\ U_a(y,\beta ,G)& =\int \frac{a(y,x)\exp ({\beta }^{T}xy)\eta (x)}{S_1(x,\beta ,G)}dx,\\ V_a(x,\beta ,G)& =\frac{1}{S_1(x,\beta ,G)}\int \frac{U_a(y,\beta ,G)}{U_1(y,\beta ,G)}\exp ({\beta }^{T}xy)dG(y),\end{array}$$

$$E\{a(Y,X)\mid X=x\}=\frac{S_a(x,{\beta }_0,G_0)}{S_1(x,{\beta }_0,G_0)},$$

$$E\{a(X,Y)\mid Y=y\}=\frac{U_a(y,{\beta }_0,G_0)}{U_1(y,{\beta }_0,G)},$$

$$E[E\{a(X,Y)\mid Y\}\mid X=x]=V_a(x,{\beta }_0,G_0).$$

$$m_a(y,x,\beta ,G)=a(y,x)-\frac{S_a(x,\beta ,G)}{S_1(x,\beta ,G)}-\frac{U_a(y,\beta ,G)}{U_1(y,\beta ,G)}+V_a(x,\beta ,G),$$

$$\begin{array}{rl}{\hat{S}}_a(x,\beta ,p)& =\sum _{k=1}^{K}a(Y_{(k)},x)\exp ({\beta }^{T}x{Y}_{(k)})p_k,\\ {\hat{U}}_a(y,\beta ,p)& =\frac{1}{n}\sum _{i=1}^n\frac{a(y,X_i)\exp ({\beta }^T{X}_{i}y)}{{\hat{S}}_1(X_i,\beta ,p)},\\ {\hat{V}}_a(x,\beta ,p)& =\frac{1}{{\hat{S}}_1(x,\beta ,p)}\sum _{k=1}^K\frac{{\hat{U}}_a(Y_{(k)},\beta ,p)}{{\hat{U}}_1(Y_{(k)},\beta ,p)}\exp ({\beta }^{T}x{Y}_{(k)})p_k.\end{array}$$

$$\frac{1}{n}\sum _{i=1}^n{m}_{n,a}\{Y_i,X_i,\beta ,\hat{p}(\beta )\}=0,$$

$$l(y,x;p,\beta )=\sum _{k=1}^K{1}_{\left\{y=Y_{(k)}\right\}}\left\{{\beta }^{T}x{Y}_{(k)}+\log (p_k)\right\}\sum _{k=1}^K{p}_k\exp \{{\beta }^{T}x{Y}_{(k)}\}.$$

$${\Psi }_n(\beta )\equiv \frac{1}{n}\sum _{i=1}^n\left[a(Y_i,X_i)-\frac{{\hat{S}}_a\{X_i,\beta ,\hat{p}(\beta )\}}{{\hat{S}}_1\{X_i,\beta ,\hat{p}(\beta )\}}-\frac{{\hat{U}}_a\{Y_i,\beta ,\hat{p}(\beta )\}}{{\hat{U}}_1\{Y_i,\beta ,\hat{p}(\beta )\}}+{\hat{V}}_a\{X_i,\beta ,\hat{p}(\beta )\}\right]=0.$$

$${\mathcal{S}}_a\equiv$$

$${\mathcal{S}}_a/{\mathcal{S}}_1\equiv$$

$$\begin{array}{rl}\underset{\beta \in \mathcal{B},p\in \mathcal{P}}{\sup }\Vert {\mathbb{P}}_n{S}_a(X,\beta ,p)-E\left\{S_a(X,\beta ,p)\right\}\Vert & \to 0,\\ \underset{\beta \in \mathcal{B},p\in \mathcal{P}}{\sup }\Vert {\mathbb{P}}_n\frac{S_a(X,\beta ,p)}{S_1(,\beta ,p)}-E\left\{\frac{S_a(X,\beta ,p)}{S_1(X,\beta ,p)}\right\}\Vert & \to 0,\end{array}$$

$$\underset{\begin{array}{c}\beta \in \mathcal{B},p\in \mathcal{P}\\ y\in \{Y_{(1)},\dots ,Y_{(k)}\}\end{array}}{\sup }\Vert {\hat{U}}_a(y,\beta ,p)-U_a(y,\beta ,p)\Vert$$

$$=\underset{\begin{array}{c}\beta \in \mathcal{B},p\in \mathcal{P}\\ y\in \{Y_{(1)},\dots ,Y_{(k)}\}\end{array}}{\sup }\Vert {\mathbb{P}}_n\frac{a(y,X)\exp ({\beta }^{T}Xy)}{S_1(X,\beta ,p)}-E\left[\frac{a(y,X)\exp ({\beta }^{T}Xy)}{S_1(X,\beta ,p)}\right]\Vert \to 0.$$

$$\begin{array}{rl}& \underset{\beta \in \mathcal{B},p\in \mathcal{P}}{\sup }\Vert {\mathbb{P}}_n\frac{{\hat{U}}_a(Y_i,\beta ,p)}{{\hat{U}}_1(Y_i,\beta ,p)}-E\left\{\frac{U_a(Y,\beta ,p)}{U_1(Y,\beta ,p)}\right\}\Vert \to 0.\end{array}$$

Full text

\startlocaldefs\setattribute

tablecaptionshape \setattributetablename skip.

T1This research was supported by Grant No. 2016126 from the United States-Israel Binational Science Foundation (BSF). We would like to thank Anastasios Tsiatis and Eric J. Tchetgen Tchetgen for very helpful discussions.

and

1 Introduction

Recently, Luo and Tsai (2012) proposed a semi-parametric proportional likelihood ratio model that extends generalized linear models. The model assumes that the joint distribution of the response $Y$ and the $q\times 1$ covariate vector $X$ is

[TABLE]

where $\beta\in\mathbb{R}^{q}$ is the Euclidean parameter of interest and $G(y)$ and $\eta(x)$ are the nuisance parameters. Here $G(y)\equiv P(Y\leq y\mid X=0)$ is a baseline distribution function with density function $g(y)$ with respect to some dominating measure $\nu$; and $\eta(x)\equiv p_{X}(x)$ is the density of $X$ with respect to some dominating measure. A comprehensive discussion of the model interpretation can be found in Luo and Tsai (2012), Chan (2013), and references therein. Semi-parametric maximum likelihood estimators of $\beta$ and $G$ were given by Luo and Tsai (2012) and the convergence of their iterative estimation algorithm was proved by Davidov and Iliopoulos (2013).

Chan (2013) showed that under a certain missing-data model (which will be discussed later), or when $Y$ is subject to doubly-random truncation, the above model is invariant with respect to $\beta$, but not necessarily with respect to $G$. By cleverly using the pairwise pseudo-likelihood of Liang and Qin (2000), Chan (2013) presented a pseudo-score equation for $\beta$ that is free of the functional parameter $G$, and such that $\beta$ is consistently estimated. This estimator is computationally efficient, but not statistically efficient. Estimation procedures based on right-censored data and on longitudinal data were proposed by Zhu (2014) and Luo (2015), respectively.

The proportional likelihood ration model (1) is a special case of the prospective likelihood discussed by Chen (2003) in the context of outcome-dependent samples. Chen (2003, Eq. 5) considered the model

[TABLE]

where

[TABLE]

and $(y_{0},x_{0})$ is a sample point. In this model, $g(y)=p_{Y\mid X,R}(y|x_{0},1)$ and $R$ is the sampling indicator, taking value $1$ if it is included in the sample and [math] otherwise. The function $\rho$ is referred to as the generalized odds ratio (Liang and Qin, 2000). For generalized linear models with canonical link function, and taking $(y_{0},x_{0})\equiv(0,0)$, $\rho$ has the from $\exp\{\beta^{T}xy\}$. Chen (2007) and Tchetgen Tchetgen et al. (2010) further extended (2) by allowing conditioning on an additional covariate vector $Z$. They study the nuisance tangent space for this model. They also considered three parametric models, for $\rho$, $p_{Y\mid X,Z}(y\mid x_{0},z)$, and $p_{X\mid Y,Z}(x\mid y_{0},z)$. They then proposed doubly robust estimators which are consistent when the model for $\rho$ is correct and either $p_{Y\mid X,Z}(y\mid x_{0},z)$ or $p_{X\mid Y,Z}(x\mid y_{0},z)$ is correctly specified.

The contribution of this work is twofold. First, it provides a comprehensive semi-parametric analysis of the Luo and Tsai’s proportional likelihood model (1) above. The proofs involve projections in Hilbert spaces and solving integral equations. Second, the semi-parametric theory we develop in this work yields an omnibus estimation procedure including problems previously studied separately, such as missing data, doubly-truncated or censored data. Under the missing-data setting, the proposed estimation approach is not limited to the specific missing-data model of Chan (2013). Moreover, in contrast to Chen (2007), Tchetgen Tchetgen et al. (2010), and Chan (2013), where only $\beta$ is estimated, we nonparametricly estimate $G$ as well, which makes our approach useful also for prediction. The utility of our novel estimation procedure is demonstrated via extensive simulation study. Efficient implementation of the proposed estimation procedure, as well as of that of Luo and Tsai (2012) and Chan (2013), is implemented in the R package PLR (which can be freely downloaded from https://github.com/yairgoldy/PLR).

2 Semi-parametric Analysis

The proportional likelihood ratio semi-parametric model can be written as the set of densities

[TABLE]

where $g(y)$ and $\eta(x)$ are the nuisance parameters. The respective true values of the parameters are denoted by $\beta_{0}$, $g_{0}(\cdot)$, and $\eta_{0}(\cdot)$. In the sequel, we assume the standard smoothness and regularity conditions (see, for example, Newey, 1990, Definition A1). Let $\mathcal{H}$ denote the tangent space for $\mathcal{P}$, where $\mathcal{H}$ is the Hilbert space of all $q$-dimensional random functions $a(Y,X)$ that satisfy $E\{a(Y,X)\}=0$ and have finite variance, equipped with the inner product $\langle a_{1},a_{2}\rangle=E(a_{1}^{T}a_{2})$. Let $\Lambda\subset\mathcal{H}$ be the nuisance tangent space with respect to the parameters $g_{0}(\cdot)$ and $\eta_{0}(\cdot)$ (see Tsiatis, 2006, Chapter 4, Defintion 1). The marginal density of $Y$ is

[TABLE]

where $g(y)$ is the density of $Y$ given $X=0$.

We start by calculating the nuisance parameters tangent space, motivated by Theorem 4.2 of Tsiatis (2006) which states that the influence function of any asymptotically linear and regular (RAL) estimator is orthogonal to the nuisance tangent space. We then show how to calculate the projection of any score function on the nuisance tangent space. As a result, we are able to provide an explicit representation of the efficient score, which is the projection of the score function with respect to $\beta$ on the orthogonal complement of the nuisance tangent space (see Tsiatis, 2006, Definition 4.2). While in practice the projections are difficult to compute since they are an infinite sequences of alternating expectations, based on an approximately-projected scores, we provide a novel family of estimators.

Lemma 1**.**

Let $\Lambda_{1}$ and $\Lambda_{2}$ be the nuisance tangent spaces with respect to $g(\cdot)$ and $\eta(\cdot)$, respectively. Then,

[TABLE]

The proof of Lemma 1 is provided in Appendix 1.1. The nuisance parameters $\eta(\cdot)$ and $g(\cdot)$ are variationally independent, that is, any choice of $\eta$ and $g$ results in a density in the model $\mathcal{P}$ (see definition at Tsiatis, 2006, page 53). Moreover, we have the following result.

Lemma 2**.**

The spaces $\Lambda_{1}$ and $\Lambda_{2}$ are orthogonal.

Proof.

For every $h(Y)-E\{h(Y)\mid X\}\in\Lambda_{1}$ and $\alpha(X)\in\Lambda_{2}$,

[TABLE]

as needed. ∎

By Theorem 5.2 of Tsiatis (2006), the projection of a function $a(Y,X)$ on $\Lambda$, denoted by $\Pi(a\mid\Lambda)$, can be written as $\Pi(a\mid\Lambda)=\Pi(a\mid\Lambda_{1})+\Pi(a\mid\Lambda_{2})$. In the following, the projection of $a$ on each nuisance tangent space is computed. Let

[TABLE]

Define the linear operator $\mathcal{B}$ by $\mathcal{B}\{f(Y)\}=E[E\{f(Y)|X\}|Y]$, and let $I:\mathcal{H}\mapsto\mathcal{H}$ denotes the identity mapping.

Theorem 3**.**

The respective projections of $a(Y,X)$ on $\Lambda_{1}$ and $\Lambda_{2}$ are

[TABLE]

See the detailed proof in the Appendix. The following statement is a direct consequence of Theorem 3, and the cornerstone for generating the proposed omnibus estimation procedure.

Corollary 4**.**

Let $\Lambda^{\bot}$ be the orthogonal complement of the nuisance tangent space $\Lambda$ in $\mathcal{H}$. Then

[TABLE]

Similar result was also obtained by Chen (2007, see the Appendix there) who used the prospective and retrospective likelihoods to shows that the space $\Lambda$ can be written as an intersection between two linear subspaces of $\mathcal{H}$. He then used von-Neumann projection theorem (Bickel et al., 1993, Theorem A.4.1) to calculate the projection. This is different from the above theorem, as here we explicitly compute the nuisance tangent spaces with respect to $g(\cdot)$ and $\eta(\cdot)$.

Let $S_{\beta}(Y,X)$ be the score function for $\beta$, namely, the derivative of $\log p_{Y,X}\{y,x;\beta,g_{0}(\cdot),\eta_{0}(\cdot)\}$ with respect to $\beta$ evaluated at the true parameter value $\beta_{0}$. The efficient score function, $S_{\mathrm{eff}}(Y,X)$, is defined as the projection of $S_{\beta}(Y,X)$ on $\Lambda^{\perp}$.

Lemma 5**.**

The efficient score $S_{\mathrm{eff}}(Y,X)$ equals

[TABLE]

See proof in the Appendix. Since the projection of $a(Y,X)$ on $\Lambda_{1}$ is an infinite series where the norm of subsequent terms decrease, we approximate the projection onto the nuisance tangent space by using only the first few terms. Our proposed approximated efficient score is defined by

[TABLE]

The asymptotic and finite-sample properties are studied in the following sections.

3 Estimation

In this section we propose a family of estimators for the parameter of interest $\beta_{0}$ using the theory developed in Section 2. An estimator $\widehat{\beta}$ for $\beta_{0}$ is called asymptotically linear if there exists a $q$-dimensional random vector $\varphi(Y,X)$, such that $E\{\varphi(Y,X)\}=0^{q\times 1}$, and $n^{1/2}(\widehat{\beta}-\beta_{0})=n^{-1/2}\sum_{i=1}^{n}\varphi(Y_{i},X_{i})+o_{p}(1)$ such that $E\{\varphi(Y,X)\varphi(Y,X)^{T}\}$ is finite and nonsingular (Tsiatis, 2006, Chapter 3). By Theorem 4.2 of Tsiatis (2006), if $\varphi$ is an influence function for an estimator $\widehat{\beta}$, then $\varphi$ is orthogonal to the nuisance tangent space $\Lambda$, or more formally,

[TABLE]

By Theorem 4.3 of Tsiatis (2006), every regular asymptotically linear (RAL) estimator for $\beta_{0}$ has a unique influence function. Therefore, we propose a family of Z-estimators, based on their influence functions using the fact that the influence functions must lie in $\Lambda^{\bot}$.

Fix any function $a(Y,X)$, $a:\mathbb{R}\times\mathbb{R}^{q}\mapsto\mathbb{R}^{q}$; then $a(Y,X)-E\{a(Y,X)\}\in\mathcal{H}$. Let

[TABLE]

The function $m_{0,a}(y,x)$ is an approximated projection of $a$ on $\Lambda^{\bot}$. Let

[TABLE]

where $S_{1}$ and $U_{1}$ are the functions $S_{a}$ and $U_{a}$, respectively, for $a\equiv 1$. Using the definition of the density $g(\cdot)$ and some algebraic manipulations, it can be shown that

[TABLE]

and

[TABLE]

Let

[TABLE]

and note that, by the above discussion, $m_{0,a}(y,x)=m_{a}(y,x,\beta_{0},G_{0})$. The function $m_{a}(y,x,\beta,G)$ is used for defining the estimating equations.

Suppose that we observe independent and identically distributed random pairs $(Y_{1},X_{1}),\ldots,(Y_{n},X_{n})$ from the distribution function $p_{Y,X}\{y,x;\beta_{0},g_{0}(\cdot),\eta_{0}(\cdot)\}$. Let $Y_{(1)},\ldots,Y_{(K)}$ be the ordered distinct observed values of $Y$. For a fixed value of $\beta$, by Theorem 2 of Luo and Tsai (2012), the profile likelihood maximizer for $G$, denote by $\widehat{G}_{\beta}$, has jumps $\widehat{p}(\beta)=\{\widehat{p}_{1}(\beta),\ldots,\widehat{p}_{K}(\beta)\}^{T}$ at $Y_{(1)},\ldots,Y_{(K)}$. For a vector of probabilities $p\in\mathbb{R}^{K}$, write

[TABLE]

Then, for every function $a(y,x)$, we propose the following estimating equation for $\beta$:

[TABLE]

where $m_{n,a}$ is defined similarly to $m_{a}$ in (4) by replacing $S_{a}$, $U_{a}$ and $V_{a}$ with $\widehat{S}_{a}$, $\widehat{U}_{a}$, and $\widehat{V}_{a}$, respectively. An estimator of $G$ is then obtained by taking $\widehat{G}=\widehat{G}_{\widehat{\beta}}$. Note that, by Lemma 5, choosing $a(Y,X)=XY$ in $m_{a}$ yields estimating equations which are based on the approximated efficient score.

4 Asymptotic Results – the discrete setting

Consider a discrete random variable $Y$ with finite support, and assume

1. (A1)

$X$ takes values in a compact set $\mathcal{X}\subset\mathbb{R}^{q}$, and $\beta_{0}$ is an interior point of a bounded set $\mathcal{B}\subset\mathbb{R}^{q}$. 2. (A2)

The function $E[m_{a}\{Y,X,\beta,p^{*}(\beta)\}]$ has a unique zero at $\beta_{0}$, and its derivative with respect to $\beta$ is invertible at $\beta_{0}$, where $p^{*}(\beta)$ is the limit of $\widehat{p}(\beta)$.

Let $p_{0}\in\mathbb{R}^{K}$ be the vector of true probabilities of $(Y_{(1)},\ldots,Y_{(K)})$ given $X=0$.

Theorem 6**.**

Under the assumptions above, $\widehat{\beta}$ and $\widehat{p}(\widehat{\beta})$ are consistent estimators for $\beta_{0}$ and $p_{0}$, respectively. Moreover, both $n^{1/2}(\widehat{\beta}-\beta_{0})$ and $n^{1/2}\{\widehat{p}(\widehat{\beta})-p_{0}\}$ converge to mean-zero Gaussian vectors.

The asymptotic variance of $\widehat{\beta}$ can be estimated empirically by standard estimating-equations tools. However, the computation is rather complex. Instead, a bootstrap approach is recommended, which is justified by Kosorok (2008, Theorem 13.4). We conjecture that Theorem 1 holds for general distributions of $Y$. Proving this is challenging. A typical first step is to show that uniformly in $\beta$, the profile likelihood maximizer $\widehat{G}_{\beta}$ converges to some limit $G_{\beta,0}$. In other words, one needs to show that uniformly in $\beta$, the random process $\operatornamewithlimits{argmax}\frac{1}{n}\sum_{i=1}^{n}\text{pl}(Y_{i},X_{i},\beta,G)$ converges to a fixed limit, where pl is the profile likelihood, and the maximization is taken over all step distribution functions with jumps at the sample points. However, the maximizer of the profile log-likelihood is given only implicitly as a solution of a nonlinear set of equations (Luo and Tsai, 2012, Theorem 2). Since no explicit solution is given for the maximizer, standard empirical process techniques are difficult to employ. This is different from proofs such as those in Murphy et al. (1997) and Luo and Tsai (2012), that use nonparametric maximum likelihood, since their proof requires convergence only at a the value of the true parameters $(\beta_{0},G_{0})$. This is also different from the locally semi-parametric proofs of Chen (2007) and Tchetgen Tchetgen et al. (2010) as a parametric model for $g$ is assumed.

Proof of Theorem 6.

Part 1: Convergence of $\widehat{p}(\beta)$ to a limit $p^{*}(\beta)$. As explained in the proof of Theorem 2 of Luo and Tsai (2012), for each fixed $\beta$, $\widehat{p}(\beta)$ is obtained by maximization the likelihood (2.2) of Vardi (1985). Note that for $\beta\neq\beta_{0}$, this maximization is carried out with respect to a misspecified model. Indeed, since $K$ is fixed, the log likelihood $l$ of one observation $(y,x)$ for a fixed $\beta$ is

[TABLE]

By Theorem 3.2 of White (1982), for every fixed $\beta$, $n^{1/2}\{\widehat{p}(\beta)-p^{*}(\beta)\}$ converges to a mean-zero Gaussian vector. Note that $l(y,x;\beta,p)$ and its first and second derivatives are all continuous function of $\beta$, and as a result of Theorem 3.2 of White (1982), the limit $p^{*}(\beta)$ is also continuous in $\beta$.

Part 2: Consistency. The outline of the consistency proof is as follows. First, we define the $\widehat{\beta}$ as a zero of an estimating equation $\Psi_{n}(\beta)=0$. We show that $\Psi_{n}(\beta)$ converges uniformly to a function $\Psi(\beta)$ which has a unique zero at $\beta_{0}$ and has the property that if $\{\beta^{(m)}\}_{m=1}^{\infty}$ is any sequence for which $\Psi\left(\beta^{(m)}\right)\rightarrow 0$, then $\beta^{(m)}\rightarrow\beta_{0}$. By Theorem 2.10 of Kosorok (2008), this proves consistency of $\widehat{\beta}$ to $\beta_{0}$.

By (6), the estimating equation is

[TABLE]

Let $\Psi(\beta)\equiv E\left[m_{a}(Y,X,\beta,p^{*}(\beta))\right]$. By Assumption (A2), $\Psi(\beta)$ has a unique zero at $\beta_{0}$. Since $p^{*}(\beta)$ is continuous and $E\{m_{a}(Y,X,\beta,p)\}$ is also continuous in both $\beta$ and $p$, so is $\Psi(\beta)$ as a composition of continuous functions. Hence, for any sequence $\{\beta^{(m)}\}$, if $\Psi(\beta^{(m)})\rightarrow 0$, then $\beta^{(m)}\rightarrow\beta_{0}$.

We now prove that $\Psi_{n}(\beta)$ converges uniformly to $\Psi(\beta)$. Define

[TABLE]

where $S_{1}$ and $S_{a}$ are defined in (5), and $\mathcal{P}\equiv\left\{p\in\mathbb{R}^{K}:p_{k}\geq 0,\,\,\sum_{k=1}^{K}p_{k}=1\right\}$. By Corollary 9.32 of Kosorok (2008), the classes $\mathcal{S}_{1}$, $\mathcal{S}_{a}$, and $\mathcal{S}_{a}/\mathcal{S}_{1}$, are Donsker since by Assumption (A1), $\beta$, $p$ and $x$ are bounded, the exponent function is Lipschitz on compact sets, and the function $a$ is bounded. Hence

[TABLE]

where $\mathbb{P}_{n}$ is the empirical measure such that for every function $f$, $\mathbb{P}_{n}f(Y,X)\equiv\sum_{i=1}^{n}f(Y_{i},X_{i})$. By Assumption (A1), $S_{1}$ is uniformly bounded from below by a positive constant. Using the same argument as above, $\mathcal{U}_{a}\equiv\large\{a(y,x)\exp(\beta^{T}xy)/S_{1}(x,\beta,p)\large\}$ is Donsker. Hence,

[TABLE]

Applying Corollary 9.32(iv) of Kosorok (2008) to the classes $\mathcal{U}_{1}$ and $\mathcal{U}_{a}$ yields that the quotient $\mathcal{U}_{a}/\mathcal{U}_{1}$ is also Donsker. Hence, one can show that

[TABLE]

Similar arguments shows that $\mathcal{V}_{a}\equiv\{V_{a}(y,\beta,p):\beta\in\mathcal{B},p\in\mathcal{P}\}$ is also Donsker and that

[TABLE]

Consequently, by the definitions of $m_{n,a}$ and $m_{a}$, and by Eqs (7), (8), and (9),

[TABLE]

which concludes the consistency proof.

Part 3: Normality. For any function $f(y,x,\beta,p)$ define $D_{\beta}f(y,x,\beta,p)\equiv\frac{\partial}{\partial\beta}f(y,x,\beta,p)$ and define $D_{p}f(y,x,\beta,p)$ similarly. Define

[TABLE]

and define $D_{0}m_{a}(y,x,\beta,p^{*}(\beta))$ similarly. By using similar arguments to (10), we get

[TABLE]

We have

[TABLE]

since $\widehat{p}(\beta_{0})-p_{0}=O_{p}(n^{-1/2})$ by Part 1 above. Write

[TABLE]

Then

[TABLE]

where

[TABLE]

and thus behaves like a V-statistic up to an $o_{p}(1)$ term. Using similar arguments, one can show that

[TABLE]

converges to a mean-zero Gaussian vector. Hence, by (11),

[TABLE]

Multiplying both sides of this equation by $-E\left[D_{0}m_{a}\left\{Y,X,\beta_{0},\widehat{p}(\beta_{0})\right\}\right]^{-1}$, and using the Donsker property for $\mathcal{S}_{a}/\mathcal{S}_{1}$, $\mathcal{U}_{a}/\mathcal{U}_{1}$ and $\mathcal{V}_{a}$, the fact that $\widehat{p}(\beta_{0})-p_{0}$ converges to a mean-zero Gaussian vector, and (12), we obtain that $n^{1/2}(\widehat{\beta}-\beta_{0})$ converges to a Gaussian random vector. Finally,

[TABLE]

which converges to a mean-zero Gaussian by Part 1 and the argument above.∎

5 Incomplete and Sampling-Biased Data

So far we assumed that the data $(Y_{1},X_{1}),\ldots,(Y_{n},X_{n})$ are fully observed and identically distributed. In the literature, the proportional likelihood model (1) with incomplete data was considered in a case by case manner. For example, Chan (2013) shows how to handle missing data when the probability of missingness has a specific form, namely,

[TABLE]

where $R$ is the indicator for non-missing data, and $h_{1}$ and $h_{2}$ are arbitrary functions. He also considers the double-truncation setting when the truncation is independent of both $Y$ and $X$. Zhu (2014) discusses the right-censored setting when the censoring variable $C$ is independent of the pair $(Y,X)$. Other settings, such as selection-biased data where the randomization is not proper, were not studied so far.

The estimating equation (6) enables us to provide an omnibus solution for all the problems that discussed aboves. Indeed, when the selection probabilities are known or can be estimated, and similarly, when the censoring or truncation probabilities can be estimated, one can use the inverse weighing methods (Robins et al., 1994). In the missing data and censoring settings, let $R_{i}$ be an indicator equals one for complete observations. Let $W_{i}=\mbox{pr}(R_{i}=1\mid Y_{i},X_{i})$, and let $\widehat{W}_{i}$ be a consistent estimator of $W_{i}$. For the sampling-biased setting, let $W_{i}$ be the sampling probability of observation $i$ and $R_{i}\equiv 1$. For a fixed $\beta$, let $\widehat{p}_{W}(\beta)$ be the weighted-profile-likelihood estimator obtained as the maximizer of

[TABLE]

Then, the solution $\widehat{\beta}_{W}$ of the estimating equation

[TABLE]

is a consistent estimator of $\beta_{0}$. Moreover, $\widehat{p}_{W}(\widehat{\beta}_{W})$ is a consistent estimator of $G$. The finite-sample performance of this estimator for the missing-data setting is demonstrated in Section 6.

6 Simulation Study

We compare our method to two existing methods: the MLE of Luo and Tsai (2012) and the pseudo-likelihood method of Chan (2013). The two scenarios of Luo and Tsai (2012) were considered. Specifically, the covariate vector consists of $X=(X_{1},X_{2})^{T}$, where $X_{2}$ follows a zero-mean normal distribution with standard deviation 0.5, and given $X_{2}$, $X_{1}$ follows the Bernoulli distribution with success probability $\exp(1-X_{2})/\{1+\exp(1-X_{2})\}$. The value of the true parameters are $\beta=(\beta_{1},\beta_{2})^{T}=(-1,-1)^{T}$. In Setting 1, $Y$ is continuous and the baseline density is defined by

[TABLE]

where $\Phi$ is the standard normal cumulative distribution function. In Setting 2, $Y$ is discrete with

[TABLE]

Each setting consists of 1000 replicates and sample sizes 100, 200, 400 and 800. We compare the bias in estimating $\beta_{1}$, $\beta_{2}$, and the distance $\int|\widehat{G}(t)-G(t)|dt$. The simulation results are summarized in Table 2 and Figure 2, in the Appendix. The proposed estimator coincides with the maximum likelihood estimator of Luo and Tsai (2012), and behaves similarly to that of Chan (2013).

Two additional settings, with missing covariates, are considered. In both settings, the full data were generated as in Setting 1. The probability of observing complete data is

[TABLE]

in Setting 3, and

[TABLE]

in Setting 4. Note that the missing probability in Setting 3 follows (13) and hence can be consistently estimated by ignoring the missing observations. This is no longer true for Setting 4. The simulation results are summarized in Table 1 and Figure 1. While all three methods work similarly in Setting 3, only the proposed method succeeds in estimating $\beta$ consistently. Moreover, the bias of the proposed method in estimating $G$ converges to zero, while for the other two methods the bias converges to a constant.

Acknowledgement

This research was supported by Grant No. 2016126 from the United States-Israel Binational Science Foundation (BSF). We would like to thank Anastasios Tsiatis and Eric J. Tchetgen Tchetgen for very helpful discussions.

Appendix

Proof of Lemma 1.

Assertion (ii) follows from  Tsiatis (2006), Theorem 4.6. For assertion (i), consider the parametric submodel

[TABLE]

where $\gamma\in\mathbb{R}^{q}$ is the nuisance parameter and $h(Y)$ is $q$-dimensional vector-valued bounded function. Clearly, the true model is obtained for $\gamma=0$. Moreover, $g_{0}(y)\exp\{h(y)^{T}\gamma\}$ is indeed a density of $Y$. The score function with respect to this submodel is given by

[TABLE]

We have demonstrated that any element in $\Lambda_{1}$ defined above is an element of a parametric submodel nuisance tangent space. Therefore, to complete the proof we need to show that the linear space spanned by the score vector with respect to $\gamma$ for any parametric submodel is contained in $\Lambda_{1}$. The log-density with respect to a parametric submodel can be written as

[TABLE]

Taking the derivative with respect to the parametric submodel $\gamma$ and substituting the true value of the parameter, denoted by $\gamma_{0}$, we obtain

[TABLE]

Multiplying the score $S_{\gamma}$ by a conformable matrix yields an element of $\Lambda_{1}$, which concludes the proof. ∎

Proof of Theorem 3.

The second assertion follows from Lemma 4.3 of Tsiatis (2006). For the first assertion, since $\Lambda_{1}\bot\Lambda_{2}$, it is enough to first project $a(Y,X)$ on $\Lambda_{2}^{\bot}$ and then on $\Lambda_{1}\subseteq\Lambda_{2}^{\bot}$. By assertion (ii), $\Pi\{a(Y,X)\mid\Lambda_{2}^{\bot}\}=a(Y,X)-E\{a(Y,X)\mid X\}$. Thus it is enough to find the projection of functions of the form $a(Y,X)-E\{a(Y,X)\mid X\}$ on $\Lambda_{1}$, for functions $a$ such that $E\{a(Y,X)\}=0$. Since all functions in $\Lambda_{1}$ are of the form $h(Y)-E\{h(Y)\mid X\}$ for some $h(Y)$, we would like to find a function $h^{*}(Y)$ such that

[TABLE]

for all $h(Y)$. Since

[TABLE]

it is enough to find $h^{*}(Y)$ such that

[TABLE]

for all $h(Y)$ in $\Lambda_{1}$. This implies that

[TABLE]

Equivalently, we would like to find $h^{*}(Y)$ that solves the integral equations

[TABLE]

where the operators $I$, $\mathcal{K}$ and $\mathcal{B}$ are defined in Section 2.

We now show that $\mathcal{B}$ is a contraction operator, that is $\left\|\mathcal{B}\{a(Y,X)\}\right\|\leq 1-\epsilon$ for some $\epsilon>0$ for all functions $a(Y,X)$ such that $\|a(Y,X)\|=1$. Denote by

[TABLE]

By Tsiatis (2006, Theorem 4.6), $\mathcal{B}\{a(Y,X)\}=\Pi\{\Pi(a(Y,X)|\Lambda_{X})|\Lambda_{Y}\}$. For any $a(Y,X)$ such that $E[a(Y,X)]=0$, and $\|a(Y,X)\|=1$, by the Pythagorean theorem,

[TABLE]

Assume that there is no positive $\varepsilon$ such that $\left\|\mathcal{B}\{a(Y,X)\}\right\|\leq(1-\epsilon)\left\|a(Y,X)\right\|$. Then, there is a sequence $\{a_{n}(Y,X)\}_{n=1}^{\infty}$ such that $\|a_{n}(Y,X)\|=1$ and $\lim_{n\rightarrow\infty}\|a_{n}(Y,X)\|=1$. By Alaoglu’s theorem (Weidmann, 2012), every bounded sequence contains a weakly convergent subsequence. Let $a_{0}$ be a limit of such a subsequence. Since $\mathcal{B}$ is a projection operator, $\|\mathcal{B}\{a_{0}(Y,X)\}\|=\|a_{0}(Y,X)\|$ which implies that $a_{0}(Y,X)$ is a function only of $Y$. Note that

[TABLE]

Hence $\|\Pi(a(Y,X)|\Lambda_{X})\|=1$ which implies that $a_{0}(Y,X)$ is a function only of $X$. Since the only function that can satisfies both conditions is a constant function, and since this function needs to have zero expectation, we arrive at a contradiction.

By Tsiatis (2006, Lemma 10.5), since $\mathcal{B}$ is a contraction operator, $h^{*}(Y)=\sum_{j=0}^{\infty}\mathcal{B}^{j}\mathcal{K}\{a(Y,X)\}$, which concludes the proof. ∎

Proof of Lemma 5.

We first compute the score function for $\beta$:

[TABLE]

The result follows from Corollary 4. ∎