Distributed Iterative Hard Thresholding for Variable Selection in Tobit Models

The analysis of left-censored data is a significant statistical focus and has attracted considerable research attention in recent years. It often arises due to the lower detection limit of an assay, posing a shared challenge across various fields such as biology, chemistry, and environmental sciences. For example, biological assays used to measure Human Immunodeficiency Virus (HIV) viral load in plasma may be limited in detecting concentrations below specific thresholds. The presence of such missingness renders commonly used linear regression methods ineffective. Furthermore, with advancements in modern data collection, these challenges often manifest in high-dimensional scenarios. In the context of HIV infection, there is a critical need to explore the association between the number of viral loads and extremely high-dimensional gene expression values. In addressing such challenges, various methodologies have been proposed in prior research. Among these, the Tobit model has proved to be valuable in modeling left-censored responses.

In recent years, significant progress has been made in the research on high dimensional censored data. A useful strategy to tackle high-dimensionality challenges involves constructing penalized estimators such as lasso-type estimators. Müller and van de Geer (2016) and Zhou and Liu (2016) provide theoretical insights into the least absolute deviation estimator with the lasso penalty. Soret et al. (2018) propose the Lasso-regularized Buckley-James least square algorithm, extending the estimator in Buckley and James (1979). Jacobson and Zou (2024) were the first to consider the high-dimensional Tobit model, which optimizes the likelihood function with a nonconvex penalty (specifically with the SCAD penalty in Fan and Li (2001)). In this paper, we concentrate on the IHT approach for variable selection. One advantage of using IHT is that the user can directly specify the number of variables to be retained, which may be useful in some scientific investigations.

IHT-style methods, which combine gradient descent with projection operations, have gained popularity in the literature for sparse recovery. Various algorithms have been proposed, such as standard IHT introduced in Blumensath and Davies (2009), GraDeS in Garg and Khandekar (2009), and Hard Thresholding Pursuit (HTP) in Foucart (2011). Jain et al. (2014) demonstrated that the IHT procedure can achieve linear convergence to the optimal solution under the conditions of strong convexity and strong smoothness in high-dimensional settings. Extending this result, Wang et al. (2023) considered nonsmooth loss functions under the less restrictive assumption of a locally positive-definite population Hessian. We will demonstrate that the linear convergence result also holds for the high-dimensional Tobit model.

Additionally, we consider a scenario where data are distributed across multiple locations and develop a distributed variant of our method. With advancements in data collection technology, gathering similar types of data from various regions has become increasingly common. Typically, directly aggregating all data at a central site faces practical hurdles related to storage, communication, and privacy concerns. To tackle these challenges, divide-and-conquer approaches are commonly employed. Many early approaches utilize one-shot methods, where estimators computed on local machines are transmitted to a central node and aggregated to form a global estimate, as in Zhang et al. (2013), Lee et al. (2017). However, such methods suffer from some drawbacks, as discussed by Jordan et al. (2019), who proposed a communication-efficient distributed algorithm. Fortunately, the gradient descent-type update in the IHT method is conducive to integration with this algorithm. Building upon this concept, we introduce a communication-efficient estimator with an IHT-type update process. To our knowledge, this represents the first exploration of the high-dimensional Tobit model in a distributed setting. Our theoretical analysis and numerical results demonstrate that this estimator’s convergence rate aligns with that achieved when pooling all data together while incurring low communication costs.

This article is structured as follows. In Section 2, we review the Tobit model and then develop both local (centralized) IHT and the distributed IHT. Section 3 is dedicated to our theoretical analysis. Section 4 draws some conclusions. All proofs are relegated to the Appendix.

We present convergence guarantees for both the local IHT and the distributed IHT in this section. In the rest of the paper, $c$ and $C$ denote positive constants independent of $n$, $d$, $M$, whose values may vary across instances. In the following, for positive sequences $a_{n}$ and $b_{n}$, $a_{n}=O\left(b_{n}\right)$ or $a_{n}\lesssim b_{n}$ means $a_{n}/b_{n}\leq c$ for some constant $c$, $a_{n}=o\left(b_{n}\right)$ means $a_{n}/b_{n}\rightarrow 0$, $a_{n}\asymp b_{n}$ means $a_{n}=O(b_{n})$ and $b_{n}=O(a_{n})$, $a_{n}\ll b_{n}$ or $b_{n}\gg a_{n}$ is the same as $a_{n}=o\left(b_{n}\right)$. We make the following assumptions.

• $\bullet$

  Assumption 1 (Truncation and signal strength): $C^{*}$ is a constant that satisfies $\left.0<C^{*}<\gamma^{*}\right.$. $\|\mbox{\boldmath$\beta$}^{*}\|_{0}=s_{0}$ and $\|\mbox{\boldmath$\beta$}^{*}\|_{2}\leq C$.

• $\bullet$

  Assumption 2 (Sparse eigenvalue condition): There exist positive constants $C_{1}$ and $\tilde{C}_{2}$ such that:

  1. (i)

     $\displaystyle C_{1}\leq E\left[\mathbf{a}^{\prime}\begin{pmatrix}I_{y>0}\mathbf{x}\mathbf{x}^{\prime}&-I_{y>0}\mathbf{x}y\\ -I_{y>0}y\mathbf{x}^{\prime}&I_{y>0}y^{2}\end{pmatrix}\mathbf{a}\right]<E\left[\mathbf{a}^{\prime}\begin{pmatrix}\mathbf{x}\mathbf{x}^{\prime}&-I_{y>0}\mathbf{x}y\\ -I_{y>0}y\mathbf{x}^{\prime}&I_{y>0}y^{2}\end{pmatrix}\mathbf{a}\right]\leq\tilde{C}_{2},$

     for all $(2s+s_{0}+1)$-sparse unit vectors $\mathbf{a}$.

  2. (ii)

     $\left(1+\sqrt{{s_{0}\over s}}\right){C_{2}-C_{1}\over C_{2}+C_{1}}<1$ holds, where $C_{2}\dot{=}\tilde{C}_{2}+{1\over\left(C^{*}\right)^{2}}$, and $s=O(s_{0})$.

• $\bullet$

  Assumption 3 (Sub-Gaussian designs): For some constant $C$, $\mathbf{x}_{i}$ are iid and $C$-sub-Gaussian, denoted as $\operatorname{subG}(C)$. That is, for all $\alpha\in\mathbb{R}^{d}$,

  $$\mathbb{E}[\exp(\alpha^{\prime}\mathbf{x})]\leq\exp\left({C^{2}\|\alpha\|_{2}^{2}\over 2}\right).$$

• $\bullet$

  Assumption 4 (Proper initialization): $\|\bm{\theta}_{int}\|_{0}\leq(s+1)$ and $\gamma_{int}\geq C^{*}$.

Remark 1: To satisfy Assumption 1, we need to use a sufficiently small $C^{*}$. The sparse eigenvalue assumption is frequently used in high-dimensional data analysis. Note that here the assumption is imposed on the population quantity so even the standard eigenvalue assumption without constraining ${\bf a}$ to be sparse is reasonable. Examining the form of the population Hessian matrix $H(\bm{\theta})$ (seen in the Appendix), it is easy to see that Assumptions 1 and 2 (i) actually imply the Hessian matrix has eigenvalues between $C_{1}$ and $C_{2}$. Assumption 3 is also a standard assumption in high-dimensional regression analysis, making it possible to derive tight bounds on different quantities (Jordan et al., 2019; Belloni and Chernozhukov, 2011). Combined with the assumption $\|\mbox{\boldmath$\beta$}^{*}\|_{2}\leq C$ and the fact that the noise is Gaussian, it also implies $y_{i}$ is sub-Gaussian. Finally, a trivial initial value satisfying assumption 4 is $\mbox{\boldmath$\theta$}_{int}=(\mathbf{0},C^{*})$.

The specific value chosen for $\eta$ is not the only option, and we opt for it for convenience in the proof. Our proof reveals that the requirement for $\eta$ and $s$ is actually that $\left(1+\sqrt{{s_{0}\over s}}\right)\max\left\{\left|1-\eta C_{1}\right|,\left|1-\eta C_{2}\right|\right\}<1$. For initialization, we recommend the Tobit-lasso estimator as it offers computational simplicity and has been shown to be consistent in Jacobson and Zou (2024).

Theorem 1 establishes the convergece of the local IHT algorithm. The following Theorem concerns the convergence of the inner loop in the distributed setting, which is very similar to Theorem 1. This is not surprising since the approximate likelihood is expected to be close to the true likelihood.

The above result shows that after the inner loop ends, we have possibly a better estimate compared to when the loop starts. Combining such bounds over multiple stages $Q$, we get the following result.

In this paper, we delved into applying the IHT algorithm to the high-dimensional Tobit model and expanded its utility to distributed environments spanning multiple machines. This marks the first effort to study the Tobit model under the distributed setting in high-dimensional scenarios. We have shown that achieving a minimax convergence rate for local and distributed estimators is feasible. Simulation results validate our theoretical findings and demonstrate the good performance of our proposed method across diverse settings. When applied to high-dimensional left-censored HIV viral load data, the IHT method also effectively produced accurate predictions.

Regarding future directions, it is possible to explore extending the IHT algorithm to decentralized distributed algorithms, such as distributed gradient descent (DGD) as introduced in Nedic and Ozdaglar (2009). DGD algorithms are extensively studied in the literature and offer advantages such as avoiding a single point of failure, addressing the limitations of centralized frameworks. Relaxing the requirement for sufficient local sample size by leveraging information from neighboring nodes in decentralized networks is also possible, as discussed in Ji et al. (2023). Exploring the combination with other distributed methods also remains an intriguing area for future research. Moreover, we can also consider a more complex class of missing data models known as Zero-Inflated models, where missing values are characterized as arising from a two-component mixture model. The presence of latent components in missing values requires a more careful analysis.