# Magnetic field driven dielectric relaxation in non-magnetic composite
medium: a low temperature study

Krishna Prasad Maity<EMAIL_ADDRESS>Department of Physics, Indian
Institute of Science, Bangalore 560012, India  Ananya Patra  Department of
Physics, Indian Institute of Science, Bangalore 560012, India  Narendra Tanty
Department of Physics, Indian Institute of Science, Bangalore 560012, India  V
Prasad  Department of Physics, Indian Institute of Science, Bangalore 560012,
India

###### Abstract

The frequency dependence of dielectric constant for composites of polyaniline
(PANI) and multi-walled carbon nanotube (MWCNT) with different degree of
functionalization is studied at low temperature (down to 4.2 K) and magnetic
field (up to 3 Tesla) applied both in parallel and perpendicular direction of
ac electric field. A relaxation phenomenon is observed in all the MWCNT/PANI
composites by applying magnetic field in both the directions, below 103 Hz.
However, PANI does not show any relaxation peak with applied magnetic field in
either direction. The relaxation peak frequency does not depend on the
strength of magnetic field but it varies with temperature and degree of
functionalization of MWCNT in composites. This relaxation phenomenon occurs
due to the inhomogeneity of the medium of two highly mismatched conductive
materials at low temperatures. We have tried to explain our results in the
light of Parish and Littlewood theory about magnetocapacitance in nonmagnetic
composite.

Université de Strasbourg, CNRS, Institut de Physique et Chimie des Matériaux
de Strasbourg, UMR 7504, 23 rue du Loess, 67000 Strasbourg, France School of
Physics and Astronomy, Tel Aviv University, 6997801 Tel Aviv, Israel

##  1 Introduction

The coupling between magnetization and dielectric constant is a very well-
known phenomenon in multiferroics which have potential application in magnetic
field sensors, multi-state memory devices 1 and so on. This coupling is
manifested by several ways in different materials. For instance,
BaTiO3-CoFe2O3 composite exhibits the magnetoelectric coupling by the strong
elastic interaction between these two phases 2. The external magnetic field
affects the magnetisation, and thus influences the dielectric properties of
these materials, known as magnetodielectric (MD) effect. However, MD effect
can be observed without the coupling between these two. The interfacial
polarization 3 at the grain boundaries can also induce the MD effect. There
are few reports where the interfacial polarization and magnetic field give the
dielectric relaxation 4. For example, La2/3Ca1/3MnO3 shows the magnetic field
dependence dielectric constant above ferromagnetic to paramagnetic transition
temperature, which is manifested by the interfacial polarization 5. Anomalous
dielectric relaxation is observed in doped silicon by applying magnetic field
perpendicular to the ac electric field 6. This effect is differentiated from
usual Debye like dielectric relaxation which occurs due to parallel magnetic
field and is explained by the effect of Lorentz force on charge polarization.
Parish and Littlewood have proposed that dielectric relaxation can be observed
in inhomogeneous medium without magnetism, where the inhomogeneity gives rise
to magnetic field dependent dielectric relaxation due to the mixing between
real and imaginary part of the dielectric constant 7. The interfacial
polarization driven MD effect in the composite of graphene with polyvinyl-
alcohol (PVA) is explained by this theory 8. The MD effect is also observed in
composites of reduced graphene oxide (RGO) with different polymers, which
occurs due to interfacial polarization, and indicates the interaction between
RGO and polymer 9.
In this article, we have studied the dielectric properties of PANI and
composite with different degree of functionalized MWCNT at low temperature
(down to 4.2 K) applying magnetic field up to 3 T. This composite is an
inhomogeneous system and becomes mixture of highly mismatched conductive
components at low temperature. The degree of functionalization of MWCNT
changes the interfacial polarization between fMWCNT and PANI. Therefore, this
system is very suitable to study the interfacial polarization in the presence
of external magnetic field. This study will serve to understand the dielectric
relaxation of non-magnetic inhomogenious medium with applied magnetic field.

##  2 Results and discussion

In Fig. 1, the variation of real part of dielectric constant
($\epsilon^{\prime}$) with frequency is shown for PANI at different
temperatures by applying magnetic field parallel to the current. We can
observe that $\epsilon^{\prime}$ decreases with increasing frequency and is
constant at high frequency. The decreasing behaviour represents the dipolar
relaxation in the system. The dipolar relaxation shift towards higher
frequency with increasing temperature. All the curves for different applied
magnetic field merge with one another at a constant temperature indicating
there is no effect of magnetic field in dielectric properties and the
electrical polarization behaviour of PANI.

Figure 1:  Variation of $\epsilon^{\prime}$ with frequency for PANI by
applying magnetic field at temperatures (a) 50 K, (b) 70 K and (c) 100 K. No
relaxation peak is observed at low frequency due to the applied magnetic
fields. [Inset graph represents the $\epsilon^{\prime}$ variation with
magnetic field, zoomed in the frequency range 40 Hz - 5 kHz at T = 50 K].

The variation of $\epsilon^{\prime}$ with frequency is shown in Fig. 2 for
different degree of functionalized fMWCNT/PANI composites by applying magnetic
field in the parallel direction of ac electric field at temperature 30 K. We
can clearly observe that a relaxation in $\epsilon^{\prime}$ appears due to
the magnetic field in all the composites at low frequency ($<$ 1000 Hz).
However, this relaxation behaviour is not present in the absence of magnetic
field (H = 0 T). This behaviour is completely dissimilar to the PANI. The
relaxation frequency does not vary with increasing magnetic field for
composites at constant temperature, however it shifts towards higher frequency
with increasing the degree of functionalization of MWNT in different
composites. The relaxation peak frequency increases from 63 Hz to 160 Hz when
the duration of functionalization increases from 0 h to 96 h. The
functionalization of MWCNT enhances the interfacial polarization and increases
the boundness of charge carrier between fMWCNT and PANI 10. The enhanced
boundness of charge carriers help to relax fast with applying magnetic field
and relaxation frequency shifts towards higher value. This coupling of
dielectric constant and magnetic field is not due to the presenece of magnetic
relaxation of either PANI or fMWCNT,since PANI manifests the paramagnetism 11
and MWCNT shows diamagnetic behaviour at low temperature 12. In our case, the
interfacial polarization between fMWCNT and PANI is affected by the magnetic
field and give rise to the dielectric relaxation at low frequency.
This magnetic field relaxation behaviour in composites can be explained by
Parish and Littlewood theory. At low temperature, PANI becomes highly
resistive but MWCNT is very conducting (conductivity difference $\approx$ 103
\- 104 S/cm). The mismatch of conductivity of two dissimilar dielectric medium
creates interfacial polarization which is affected by the magnetic field. PANI
contains relatively high and low conductive regions depending on the doping.
The charge transport between these high and low conductive regions creates the
electrical polarization and shows the finite dielectric constant. At low
temperature, PANI behaves like a homogeneous system of almost similar
conductive regions. In this case, external magnetic field does not affect the
dielectric constant and no relaxation is observed at low frequency. At high
temperature the electrode polarization might contibute to the dielectric
relaxation phenomenon and induce magnetic field dependence. For proper
comparison, we have used same silver plate electrodes for pure PANI and the
composites. However, magnetic field dependence is not observed in pure PANI.
Hence, contribution from electrode towards magnetic field dependence
relaxation can be ruled out.

Figure 2: Variation of real part of dielectric constant with frequency
applying magnetic field at temperature 30 K for (a) MWCNT/PANI (0 h), (b) 48 h
and (c) 96 h fMWCNT composite with PANI. The relaxation shifts towards higher
frequency with increasing degree of functionalization of MWCNT. [The legend of
graph (a) is also shared by (b) and (c).]

The variation of $\epsilon^{\prime}$ with frequency of 96 h fMWCNT/PANI
composite is shown in Fig. 3 for different temperatures. The relaxation
behaviour is observed at all the temperatures and relaxation frequency shifts
towards lower frequency with increasing temperature. The thermal energy of
charge carrier increases when temperature increases and deflect longer
distance due to enhanced Lorentz force, and charge carriers take longer time
to relax with magnetic field i.e. the relaxaton process become slow.

Figure 3:  Variation of real part of dielectric constant with frequency for 96
h fMWCNT/PANI composite applying magnetic field at different temperature (a)
4.2 K, (b) 10 K, (c) 20 K and (d) 30 K. The relaxation shifts towards lower
frequency with increasing temperature.

In Fig. 4 the relaxation behaviour of $\epsilon^{\prime}$ is shown at 30 K by
varying the magnetic field. The difference between dielectric constant with
magnetic field and without magnetic field, for both real and imaginary part
[$\Delta\epsilon^{\prime}=\epsilon^{\prime}(H)-\epsilon^{\prime}(0)$
similarly, $\Delta\epsilon^{\prime\prime}$] are plotted in Fig. 4 (a) and (b),
respectively. It is important to note that the relaxation frequency is not
shifted with increasing magnetic field, but the width and height of the peak
are increased. It can be speculated that the magnetic field excites more
number of charge carriers to relax at the same frequency, and does not change
the relaxation time of the charge carrier.

Figure 4:  The difference between with and without applying magnetic field of
(a) real and (b) imaginary parts of dielectric constant as a function of
frequency for 96 h fMWCNT/PANI composites at 30 K.

We have measured the dielectric relaxation behaviour of PANI and fMWCNT/PANI
composites by applying magnetic field in the perpendicular direction of
applied ac electric field and plotted in Fig. 5. No relaxation peak is
observed for PANI at low frequency for applying magnetic field up to 5 T. The
MWCNT/PANI exhibits relaxation and the peak intensity increases with
increasing magnetic field. However, the relaxation peak frequency does not
shift position with magnetic field similar to the applied field in parallel
direction. The relaxation behaviour of MWCNT/PANI composite in the both
parallel and perpendicular direction imply the inhomogeneity of the composite
system is distributed isotropically in all directions e.g. 3D nature of the
system. The MWCNT makes network in the PANI and creates interfacial
polarization which are isotropically distributed in 3D fashion. We can also
observe that another relaxation appears at higher frequency for perpendicular
field. This indicates the higher energetic relaxation in the system.

Our result is very dissimilar to the doped silicon and RGO composites with
different polymers. In doped silicon, Brook et al. observed the anomalous
dielectric relaxtion for perpendicular direction however it is usual Debye
like for parallel field 6. They have explained that polarized charge carriers
feel the Lorentz force due to applied magnetic field in perpendicular
direction and modify the polarization. The relaxation peak frequency shifts
with increasing magnetic field for a constant temperature. In our case, there
is no shift in relaxation peak by increasing magnetic field, and can’t be
explained completely by considering this model. On the other hand, Parish and
Littlewood (PL) proposed that there is possibility to observe the dielectric
relaxation [we have mensioned here ‘relaxation’ as it takes place at low
frequency though PL used the term ‘resonance’] in 2D inhomogeneous medium and
interface between different conductive materials 13, 7 by applying magnetic
field without magnetism present in the system. The dielectric relaxation is
obeserved due to the mixing of real and imaginary part of the dielectric
response occured by inhomogeneity of the medium. The relaxation peak is
observed for the condition $\beta\omega\tau$ = 1, where $\beta$ = $\mu$H
($\mu$ is the mobility of charge carrier and H is magnetic field) and $\tau$ =
$\rho\epsilon$ ($\rho$ is resistivity). Considering the resistivity variation
with temperature is like $\rho$ = $\rho_{0}$exp($\Delta$/kBT); $\Delta$ is
activation energy, the relaxation condition modifies to
ln($\frac{1}{\beta\omega\epsilon\rho_{0}}$)= $\dfrac{\Delta}{k_{B}T}$.
Assuming no variation of dielectric constant with temperature, modified
relaxation condition predicts that the relaxation peak shifts towards higher
frequency with increasing temperature for a constant magnetic field, and it
should shift towards lower frequency with increasing magnetic field for a
constant temperature.

Figure 5:  Variation of real part of dielectric constant with frequency for
(a) PANI and (b) MWCNT/PANI applying magnetic field perpendicular to the ac
electric field.

In case of our sample, the dielectric resonance peak shifts towards lower
frequency with increasing temperature for a constant magnetic field which is
the opposite behaviour of the prediction. PL model also tells that the
dielectric resonance enhances for the current flow in the perpendicular
direction of the magnetic field. In the MWCNT/PANI composites, MWCNT is
distributed in all the direction (3D nature) and finds the perpendicular
current flow at the interface with respect to the magnetic field for which we
have observed the dielectric relaxation for both parallel and perpendicular
direction.
In summary, we studied the dielectric properties of non-magnetic fMWCNT/PANI
composite with varying functionalization, temperature and magnetic field. We
have found that the magnetic field has significant effect on the origin of
dielectric relaxation. The dielectric relaxation is observed for all the
composite at low frequency in the presence of magnetic field. We have tried to
explain our results with the essence of 2D PL theory. However, the results are
not fully explained by this theory. An extended version of this theory in
3-dimension is required to understand the magnetic field dependence dielectric
relaxation. Interestingly, the dielectric relaxation due to inhomogeneity of
the system is studied extensively in this fMWCNT/PANI system. We believe this
study will help to pave the way of this phenomena in other similar composite
systems.

##  3 Experimental

The chemical functionalization of MWCNT is done by using concentrated H2SO4
and HNO3 (vol. ratio 3:1) acid mixtures. The degree of functionalization of
MWCNT is varied by changing the duration (hour; denoted as ‘h’; we have
prepared particularly 6 h, 48 h and 96 h) of functionalization 14. The
composites of MWCNT and PANI are prepared by using in-situ polymerization
technique which is explained in our previous reports 15. The samples were made
in the form of pallet with radius 10 mm and thickness 1.35 mm. The dielectric
measurements are performed keeping the samples in between two silver plates,
and put inside the Janis cryostat equipped with superconducting magnet. The
magnetic field is applied both in perpendicular and parallel directions of
applied ac electric field by altering the orientation of the sample with
respect to the external magnetic field. The dielectric measurements are
performed by using Agilent 4294A high precision impedance analyzer.
We have also characterized the samples using SEM, Raman spectroscopy, XPS in
our previous reports. It is observed that PANI coated on the surface of MWCNT
and distributed randomly in the system. The different functional groups
attached on the surface of MWCNT enhances the polar interaction between fMWCNT
and PANI.

Authors thank to Seena Mathew and Husna Jan for helping in measurements. KPM
thanks to IISc for reasearch associate fellowship.

## References

  * Scott 2007 Scott, J. F. _Science_ 2007, _315_ , 954–959
  * Zheng et al. 2004 Zheng, H. et al. _Science_ 2004, _303_ , 661–663
  * Catalan 2006 Catalan, G. _Applied Physics Letters_ 2006, _88_ , 102902
  * Maglione 2008 Maglione, M. _Journal of Physics: Condensed Matter_ 2008, _20_ , 322202
  * Rivas et al. 2006 Rivas, J.; Mira, J.; Rivas-Murias, B.; Fondado, A.; Dec, J.; Kleemann, W.; Senaris-Rodriguez, M. A. _Applied Physics Letters_ 2006, _88_ , 242906
  * Brooks et al. 2008 Brooks, J. S.; Vasic, R.; Kismarahardja, A.; Steven, E.; Tokumoto, T.; Schlottmann, P.; Kelly, S. _Phys. Rev. B_ 2008, _78_ , 045205
  * Parish and Littlewood 2008 Parish, M. M.; Littlewood, P. B. _Phys. Rev. Lett._ 2008, _101_ , 166602
  * Mitra et al. 2011 Mitra, S.; Mondal, O.; Saha, D. R.; Datta, A.; Banerjee, S.; Chakravorty, D. _The Journal of Physical Chemistry C_ 2011, _115_ , 14285–14289
  * Pradhan et al. 2020 Pradhan, S. S.; Ghosh, T. N.; Marik, A.; Raul, K. K.; Sarkar, S. K. _Bulletin of Materials Science_ 2020, _43_ , 208
  * Maity et al. 2020 Maity, K. P.; Patra, A.; Prasad, V. _Journal of Physics D: Applied Physics_ 2020,
  * Novak et al. 2010 Novak, M.; Kokanovic, I.; Babic, D.; Bacani, M. _Journal of Non-Crystalline Solids_ 2010, _356_ , 1725 – 1729
  * Ellis and Ingham 2006 Ellis, A.; Ingham, B. _Journal of Magnetism and Magnetic Materials_ 2006, _302_ , 378 – 381
  * Parish 2014 Parish, M. M. _Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences_ 2014, _372_ , 20120452
  * Maity and Prasad 2019 Maity, K. P.; Prasad, V. _Materials Research Express_ 2019, _6_ , 0850a2
  * Maity et al. 2018 Maity, K. P.; Tanty, N.; Patra, A.; Prasad, V. _Materials Research Express_ 2018, _5_ , 035034


# Nonparametric augmented probability weighting with sparsity

Xin He, Xiaojun Mao, and Zhonglei Wang X. He is with the School of Statistics
and Management, Shanghai University of Finance and Economics, Shanghai 200433,
China (e-mail:he.xin17@mail.shufe.edu.cn); X. Mao is with the School of
Mathematical Sciences, Shanghai Jiao Tong University, Shanghai 200240, China
(e-mail: maoxj@sjtu.edu.cn). Z. Wang is with Wang Yanan Institute for Studies
in Economics and School of Economics, Xiamen University, Xiamen, Fujian
361005, China (e-mail: wangzl@xmu.edu.cn). X. He and X. Mao contributed
equally to this work. Corresponding authors: Z. Wang.

###### Abstract

Nonresponse frequently arises in practice, and simply ignoring it may lead to
erroneous inference. Besides, the number of collected covariates may increase
as the sample size in modern statistics, so parametric imputation or
propensity score weighting usually leads to inefficiency without consideration
of sparsity. In this paper, we propose a nonparametric imputation method with
sparse learning by employing an efficient kernel-based learning gradient
algorithm to identify truly informative covariates. Moreover, an augmented
probability weighting framework is adopted to improve the estimation
efficiency of the nonparametric imputation method and establish the limiting
distribution of the corresponding estimator under regularity assumptions. The
performance of the proposed method is also supported by several simulated
examples and one real-life analysis.

_Keywords:_ Central limit theorem, Reproducing kernel Hilbert space,
Nonresponse, Sparse learning.

##  1 Introduction

Nonresponse is a common problem in social science and other related fields,
and simply ignoring it may lead to inefficiency or even erroneous inference
due to confounding covariates (Rosenbaum and Rubin, 1983; Qu et al., 2010;
Abadie and Imbens, 2016; Lin et al., 2018). Moreover, the number of collected
covariates is relatively large in modern statistics, which makes learning
nonresponse even more challenging (Yang et al., 2020). How to deal with the
nonresponse under the high-dimensional setup is still an open question.

Sparse learning bridges the gap between the high-dimensional data analysis and
nonresponse. It is generally believed that among the numerous covariates, only
a few of them contribute to the response of interest, known as truly
informative ones, while others are noise. Thus, a variety of sparse learning
methods have been proposed to identify those truly informative covariates
under regularity assumptions. The linear response model assumption is
popularly imposed, and various attempts have been made by designing sparse-
induced regularization (Tibshirani, 1996; Fan and Li, 2001; Zou, 2006; Shen et
al., 2012, 2013), evaluating the marginal dependence (Fan and Lv, 2008; Wang
and Leng, 2016), or checking variable robustness against added noise (Barber
and Candes, 2019). Extended methods have also been developed for nonparametric
models (Lin and Zhang, 2006; Huang et al., 2010; Fan et al., 2011). However,
all these methods require explicit model assumptions that are difficult to
validate in practice or suffer from heavy computational burden. To circumvent
this difficulty, Belloni et al. (2013) proposed a feasible lasso method, which
is similar to the adaptive lasso (Zou, 2006), for variable selection in a
partial linear model. By using machine learning algorithms to handle high-
dimensional nuisance parameters, Chernozhukov et al. (2018) proposed a
double/debiased machine learning procedure to achieve parametric convergence
rate for a low dimensional parameter. A valid double robust estimator using
lasso-type penalty is discussed by Tan (2020). Recently, kernel-based sparse
learning methods have been inspired by the fact that the gradient functions
provide an appropriate criterion to identify a general dependence structure in
a model-free fashion. Specifically, Rosasco et al. (2013) proposed a novel
learning-gradient method, which adds an empirical functional penalty on the
gradients to a standard kernel ridge regression in a reproducing kernel
Hilbert space (RKHS). Besides, Yang et al. (2016) employed pair-wise learning
to estimate the gradient functions and considered a functional group lasso
penalty to induce sparsity. It is worth pointing out that the lack of
selection consistency (Rosasco et al., 2013) and the high computational cost
(Yang et al., 2016) remain unsolved. To alleviate the difficulties, He et al.
(2021) proposed a two-step sparse learning framework, which is computationally
efficient in the sense that it only requires to fit the standard kernel ridge
regression and the selection consistency is established under regularity
assumptions. The method proposed by He et al. (2021) can be regarded as a
nonparametric joint screening approach and achieves methodological
flexibility, numerical efficiency and asymptotic consistency simultaneously.

Propensity score weighting is commonly used to handle nonresponse (Robins et
al., 1994; Wooldridge, 2007; Tan, 2010; Graham et al., 2012; Zhao et al.,
2017), but conventional methods using all covariates may lead to numerical
failure, including the lack of convergence and inefficiency, due to
overfitting. Thus, sparse assumption is often imposed to estimate the
propensity scores more efficiently (Shevade and Keerthi, 2003; Genkin et al.,
2007). A Bayesian variable selection method has been proposed by Chen et al.
(1999) for logistic regression; also see Wainwright et al. (2007), Banerjee et
al. (2008) and Ravikumar et al. (2010) for details about penalized logistic
regression models. The group lasso (Yuan and Lin, 2006) was generalized to
logistic regression model by Kim et al. (2006), and Meier et al. (2008)
proposed a more efficient group lasso algorithm than that of Kim et al.
(2006). Besides, Meier et al. (2008) also established the asymptotic
consistency of the corresponding estimator. Ning et al. (2020) proposed a
high-dimensional covariate balancing propensity score estimator, and they
validated that their proposed estimator is of parametric convergence rate and
is asymptotically normal distributed. See Tang et al. (2014) and Bertsimas and
King (2017) for a review of identifying informative covariates for logistic
regression models.

In this paper, we propose a nonparametric Augmented Inverse Probability
Weighting (AIPW) framework (Robins et al., 1994) to handle the nonresponse
under the assumption of sparsity. Inspired by the key observation that the
gradient functions provide an appropriate information of the truly informative
covariates in a model-free fashion, we employ a kernel-based sparse learning
algorithm to efficiently impute the nonresponse. It only requires to fit the
standard kernel ridge regression, which has an analytical solution, and the
gradient functions can be directly computed by the derivative reproducing
property. More importantly, the truly informative covariates can be exactly
recovered with high probability. Even though the nonparametric imputation with
sparse learning achieves consistency, its convergence rate is at most
$O_{\mathbb{P}}(n^{-1/6}\log(n))$ under regularity assumptions, so it is hard
to construct an interval estimator. To alleviate this difficulty, an AIPW
framework is adopted to improve the convergence rate of the corresponding
estimator, and a central limit theorem can be established. To achieve this
goal, certain propensity score methods for analyzing sparse data suffice under
regularity conditions; see Section 2.2 for details. The corresponding variance
estimator is also discussed. The superior performance of the proposed
nonparametric AIPW framework is also supported by the numerical comparisons
against some state-of-the-art methods in several simulated examples and one
real-life analysis.

The rest of this paper is organized as follows. Section 2 provides the
background and introduces the proposed nonparametric AIPW framework. The
theoretical properties of the corresponding estimator are established under
regularity assumptions in Section 3. Section 4 reports the numerical
experiments on the simulated and real-life examples. A brief summary is
provided in Section 5.

##  2 Method

Consider

$\displaystyle Y=f^{*}(\bm{x})+\epsilon,$ (1)

where $f^{*}(\bm{x})=\mathbb{E}(Y\mid\bm{X}=\bm{x})$ is a continuous function
of a covariate vector $\bm{x}=(x_{1},...,x_{p})^{\top}$ taking values from a
$p$ dimensional separable and compact metric space ${\cal
X}\subset\mathbb{R}^{p}$, and $\epsilon$ denotes a random noise with
conditional mean zero and bounded variance. We are interested in inferring
$\theta^{*}=\mathbb{E}(Y)$ from a random sample
$\\{(\bm{x}_{i},y_{i}):i=1,\ldots,n\\}$ generated by (1).

If $y_{1},\ldots,y_{n}$ were fully observed, the sample mean
$\widehat{\theta}_{n}=n^{-1}\sum_{i=1}^{n}y_{i}$ would be an efficient
estimator of $\theta^{*}$. However, it is generally not the case in practice,
and the response of interest suffers from nonresponse. For $i=1,\ldots,n$,
denote $\delta_{i}$ to be the response indicator of $y_{i}$, where
$\delta_{i}=1$ if $y_{i}$ is observed and 0 otherwise. For simplicity, we
assume missing at random (Rubin, 1976) for the response mechanism,

$\Pr(\delta_{i}=1\mid\bm{x}_{i},y_{i})=\Pr(\delta_{i}=1\mid\bm{x}_{i}),$ (2)

and denote $\pi^{*}(\bm{x})=\Pr(\delta=1\mid\bm{x})$.

If consistent estimators $\widehat{f}_{0}(\bm{x})$ and $\widehat{\pi}(\bm{x})$
for $f^{*}(\bm{x})$ and $\pi^{*}(\bm{x})$ are available, then an AIPW
estimator,

$\widehat{\theta}_{AIPW}=\frac{1}{n}\sum_{i=1}^{n}\left[\widehat{f}_{0}(\bm{x}_{i})+\frac{\delta_{i}}{\widehat{\pi}(\bm{x}_{i})}\\{y_{i}-\widehat{f}_{0}(\bm{x}_{i})\\}\right],$
(3)

can be applied to estimate $\theta^{*}$. More rigorously, the estimator (3) is
not an AIPW estimator unless we replace $n^{-1}$ by
$(\sum_{i=1}^{n}\delta_{i}\widehat{\pi}_{i}^{-1})^{-1}$. However, if the
response model is correctly specified, we can show that
$n^{-1}\sum_{i=1}^{n}\delta_{i}\widehat{\pi}_{i}^{-1}\to 1$ in probability
under regularity conditions. When $p$ is small and $f^{*}(\bm{x})$ is a
parametric model, standard statistical methods can be used to obtain
$\widehat{f}_{0}(\bm{x})$ and $\widehat{\pi}(\bm{x})$. As $p$ increases,
however, it is not reasonable to include all covariates to estimate
$f^{*}(\bm{x})$ and $\pi^{*}(\bm{x})$ due to the curse of dimensionality.
Moreover, the model misspecification for $f^{*}(\bm{x})$ leads to erroneous
inference.

###  2.1 Estimation of ${f}^{*}(x)$ via nonparametric sparse learning

To overcome those difficulties, we employ an efficient kernel-based sparse
learning algorithm (He et al., 2021) to estimate $f^{*}(\bm{x})$ in (1).
Denote ${\cal H}_{K}$ to be an RKHS induced by a pre-specified kernel function
$K(\cdot,\cdot)$, where $K(\cdot,\cdot):{\cal X}\times{\cal
X}\rightarrow{\mathbb{R}}$ is bounded, symmetric and positive semi-definite.
It can be shown that ${\cal H}_{K}$ associated with the kernel
$K(\cdot,\cdot)$ is the completion of the linear space spanned by
$\\{K_{\bm{x}}(\cdot):\bm{x}{\in}{\cal X}\\}$ with an inner product $\langle
K_{\bm{x}},K_{\bm{u}}\rangle_{K}=K(\bm{x},\bm{u})$ for $\bm{x},\bm{u}\in{\cal
X}$, where $K_{\bm{x}}(\cdot)=K(\bm{x},\cdot)$. Thus, ${\cal H}_{K}$ is
uniquely determined by a kernel function $K(\cdot,\cdot)$ and the reproducing
property, $\langle f,{K}_{\bm{x}}\rangle_{K}=f(\bm{x})$ for $f\in{\cal H}_{K}$
and $\bm{x}\in{\cal X}$. It is noteworthy that the RKHS induced by some
universal kernel, such as the Gaussian kernel, is fairly large in the sense
that any continuous function can be well approximated (Steinwart, 2005). Thus,
by (2), we assume

$\displaystyle f^{*}(\bm{x})=\mathop{\rm argmin}_{f\in{\cal
H}_{K}}\mathbb{E}[\delta\\{y-f({\bm{x}})\\}]^{2}.$ (4)

By (2), (4) implies that $f^{*}(\bm{x})=\mathbb{E}(Y\mid\bm{X}=\bm{x})\in{\cal
H}_{K}$.

A covariate $x_{l}$ is non-informative for $f^{*}(\bm{x})$, if and only if
$g^{*}_{l}(\bm{x})=0$ for ${\bm{x}}\in{\cal X}$ almost surely, where $x_{l}$
is the $l$th component of $\bm{x}$ for $l=1,\ldots,p$, and
$g_{l}^{*}(\bm{x})=\partial f^{*}(\bm{x})/\partial x_{l}$. Thus, the
usefulness of $x_{l}$ for estimating $f^{*}(\bm{x})$ can be measured via the
${\cal L}^{2}$-norm of $g^{*}_{l}(\bm{x})$, $\|g^{*}_{l}\|^{2}_{2}=\int_{\cal
X}\\{g_{l}^{*}(\bm{x})\\}^{2}d\rho(\bm{x}),$ where $\rho(\bm{x})$ is the
marginal distribution of $\bm{X}$. Denote ${\cal
A}^{*}=\big{\\{}l:\|g^{*}_{l}\|^{2}_{2}>0\big{\\}}$ to be the active set
containing all informative covariates associated with $f^{*}(\bm{x})$. To
estimate $g^{*}_{l}(\bm{x})$ efficiently, we consider the following derivative
reproducing property (Zhou, 2007):

$\displaystyle g^{*}_{l}(\bm{x})=\langle
f^{*},{\partial_{l}{K}_{\bm{x}}}\rangle_{K},$ (5)

where $\partial_{l}{K}_{\bm{x}}(\cdot)=\partial K(\bm{x},\cdot)/\partial
x_{l}$. Specifically, by (4)–(5), once an initial estimator of $f^{*}(\bm{x})$
is available, say $\widehat{f}(\bm{x})$, its gradient function
$g^{*}_{l}(\bm{x})$ can be estimated by (5). Without loss of generality, we
assume that the first $m$ samples are fully observed, and the subsequent
analysis is conditional on the realised sample, where
$m=\sum_{i=1}^{n}\delta_{i}$. Under (2) and the high-dimensional setup, to
obtain an initial estimator $\widehat{f}(\bm{x})$, we employ the standard
kernel ridge regression by solving

$\widehat{f}=\underset{f\in{\cal
H}_{K}}{\arg\min}\Big{[}\frac{1}{m}\sum_{i=1}^{m}\left\\{y_{i}-f(\bm{x}_{i})\right\\}^{2}+\lambda\left\|f\right\|_{K}^{2}\Big{]},$
(6)

where $\lambda>0$ is a tuning parameter controlling the model complexity and
typically goes to 0 as $m$ goes to infinity, and
$\left\|\cdot\right\|_{K}^{2}$ denotes the RKHS-norm induced by the inner
product $\langle\cdot,\cdot\rangle_{K}$; see Section 3 for details. By the
representer theorem (Mercer, 1909), the solution of (6) is of the form

$\displaystyle\widehat{f}(\bm{x})=\sum_{i=1}^{m}\widehat{\alpha}_{i}K({\bm{x}}_{i},\bm{x})=\widehat{\bm{\alpha}}^{\top}{\bm{K}}_{m}(\bm{x}),$
(7)

where
${\bm{K}}_{m}(\bm{x})=(K({\bm{x}}_{1},\bm{x}),...,K({\bm{x}}_{m},\bm{x}))^{\top}$,
and
$\widehat{\bm{\alpha}}=(\widehat{\alpha}_{1},...,\widehat{\alpha}_{m})^{\top}\subset\mathbb{R}^{m}$
are the estimated representer coefficients. That is, the representer theorem
converts the original optimization task (6) in an infinite functional space
${\cal H}_{K}$ into a $m$-dimensional vector space. By (7), the optimization
task (6) is equivalent to

$\widehat{\bm{\alpha}}=\underset{\bm{\alpha}}{\mathop{\rm
argmin}}\Big{[}\frac{1}{m}\sum_{i=1}^{m}\left\\{y_{i}-\bm{\alpha}^{\top}{\bm{K}}_{m}(\bm{x}_{i})\right\\}^{2}+\lambda\bm{\alpha}^{T}\bm{K}\bm{\alpha}\Big{]},$

and its solution is
$\widehat{\bm{\alpha}}=(\bm{K}+\lambda\bm{I})^{-1}{\bm{y}},$ where $\bm{K}$ is
an $m\times m$ matrix with $(i,j)$th entry being $K(\bm{x}_{i},\bm{x}_{j})$
and $\bm{y}=(y_{1},\ldots,y_{m})^{\top}\in\mathbb{R}^{m}$.

Once $\widehat{\bm{\alpha}}$ is obtained, the gradient function in (5) can be
estimated by

$\widehat{g}_{l}(\bm{x})=\frac{\partial\widehat{f}({\bm{x}})}{\partial
x_{l}}=\widehat{\bm{\alpha}}^{\top}{\partial_{l}{\bm{K}}_{m}({\bm{x}})},\quad(l=1,...,p),$

where
${\partial_{l}{\bm{K}}_{m}({\bm{x}})}=(\partial_{l}{K}_{\bm{x}_{1}}({\bm{x}}),\ldots,\partial_{l}{K}_{\bm{x}_{m}}({\bm{x}}))^{\top}$.
Since the marginal distribution $\rho(\bm{x})$ is seldom available, instead of
the ${\cal L}^{2}$-norm, the empirical norm $\left\|\cdot\right\|_{m}$ is
considered:

$\|\widehat{g}_{l}\|^{2}_{m}=\frac{1}{m}\sum_{i=1}^{m}\big{\\{}\widehat{g}_{l}(\bm{x}_{i})\big{\\}}^{2}=\frac{1}{m}\sum_{i=1}^{m}\big{\\{}\widehat{\bm{\alpha}}^{\top}{\partial_{l}{\bm{K}}_{m}({\bm{x}}_{i})}\big{\\}}^{2},$

and the estimated active set is $\widehat{\cal
A}_{v_{m}}=\Big{\\{}l:\|\widehat{g}_{l}\|^{2}_{m}>v_{m}\Big{\\}},$ where
$v_{m}$ is a thresholding value determined through a stability-based selection
criterion (Sun et al., 2013). Finally, we refit (6) with the selected
covariates in $\widehat{\cal A}_{v_{n}}$ to obtain the nonparametric estimator
$\widehat{f}_{0}(\bm{x})$.

It is worthy pointing out that the employed sparse learning algorithm was
originally proposed by He et al. (2021), and they only focused on the purpose
of variable selection and established the selection consistency without
considering nonresponse. Yet, we generalized their method to handle incomplete
samples in this paper and treat it as an valid intermediate estimator of our
proposed estimator. More importantly, we further established a central limit
theorem for the proposed nonparamtric estimator, which is rare and attractive
in machine learning, and a variance estimation is also provided as well; see
Section 3 for details.

###  2.2 Estimation of $\pi^{*}(x)$

There exist works to estimate $\pi^{*}(\bm{x})$ under the assumption of
sparsity, and we consider the group lasso (Meier et al., 2008) as an example
by assuming

$\mbox{logit}\\{\pi^{*}(\bm{x}_{i})\\}=\beta^{*}_{0}+\bm{x}_{i}^{\top}\bm{\beta}^{*}_{1},$

where $\mbox{logit}(z)=\log(z)-\log(1-z)$ for $z\in(0,1)$. In addition, assume
that the covariate vector can be rewritten as
$\bm{x}=(\bm{x}_{1}^{\top},\ldots,\bm{x}_{G}^{\top})^{\top}$, where
$\bm{x}_{g}\in\mathbb{R}^{\mbox{\scriptsize df}_{g}}$ contains the covariates
of the $g$th group for $g=1,\ldots,G$, and $\mbox{df}_{g}$ is the
corresponding degrees of freedom. For example, $\mbox{df}_{g}=3$ if
$\bm{x}_{g}$ corresponds to a categorical covariate with four levels, and
$\mbox{df}_{g}=1$ if $x_{g}$ is continuous; see Meier et al. (2008) for
details.

The log-likelihood estimator with a group lasso penalty is obtained by solving

$\displaystyle\widehat{\bm{\beta}}_{\lambda_{2}}$ $\displaystyle=$
$\displaystyle\mathop{\rm
argmin}_{\bm{\beta}}\Big{\\{}-l(\bm{\beta})+\lambda_{2}p(\bm{\beta})\Big{\\}},$
(8)

where
$l(\bm{\beta})=\sum_{i=1}^{n}[\delta_{i}\log\\{\pi(\bm{x}_{i})\\}+(1-\delta_{i})\log\\{1-\pi(\bm{x}_{i})\\}]$
is the log-likelihood of the response indicators,
$\bm{\beta}^{\top}=(\beta_{0},\bm{\beta}_{1}^{\top})\in\mathbb{R}^{p+1}$ with
$\beta_{0}\in\mathbb{R}$ and $\bm{\beta}_{1}\in\mathbb{R}^{p}$,
$p(\bm{\beta})=\sum_{g=1}^{G}\mbox{df}_{g}^{1/2}\lVert\bm{\beta}_{g}\rVert_{2}$
is the group lasso penalty, $\lVert\cdot\rVert_{2}$ is the Euclidean norm, and
$\bm{\beta}_{g}$ corresponds to $\bm{x}_{i,g}$ for $g=1,\ldots,G$. The block
co-ordinate gradient descent algorithm is used to obtain
$\widehat{\bm{\beta}}_{\lambda_{2}}$ in (8), and the detailed algorithm is
adjourned to Appendix A.

######  Remark 1.

Since the estimated response probability is used to improve the convergence
rate of the estimator in (3), the response model is assumed to be correctly
specified; see Qin et al. (2017) for a similar assumption. In addition to the
group lasso method (Meier et al., 2008), other penalized logistic regression
estimators (Fan et al., 2014; Ning et al., 2020) can be used to estimate the
response probability. However, to guarantee the asymptotic central limit
theorem in Theorem 2, the estimated response probability by other methods
should satisfy Lemma 3; see Section 3 for details.

##  3 Theoretical Properties

In this section, we investigate the asymptotic consistency of
$\widehat{f}_{0}(\bm{x})$ and establish the central limit theorem for the AIPW
estimator in (3) under regularity assumptions.

Denote an integral operator $L_{K}:{\cal L}^{2}({\cal
X},{\rho})\rightarrow{\cal L}^{2}({\cal X},{\rho})$ as
$L_{K}(f)(\bm{x})=\int_{\cal X}K(\bm{x},\bm{u})f(\bm{u})d\rho(\bm{u}),$ for
$f\in{\cal L}^{2}({\cal X},{\rho})$, where ${\cal L}^{2}({\cal
X},{\rho})=\\{f:\int f^{2}(\bm{x})d\rho(\bm{x})<\infty\\}$. If the RKHS ${\cal
H}_{K}$ is separable, then by the spectral theorem (Fischer and Steinwart,
2020), we have $L_{K}f=\sum_{j\geq 1}\mu_{j}\langle f,e_{j}\rangle_{2}e_{j},$
where $\\{e_{j}:j=1,2,\ldots\\}$ form an orthonormal basis of ${\cal
L}^{2}({\cal X},\rho)$, $\\{\mu_{j}:j=1,2,\ldots\\}$ are the corresponding
eigenvalues with respect to $L_{K}$, and $\langle f,g\rangle_{2}=\int_{\cal
X}f(\bm{x})g(\bm{x})d\rho(\bm{x})$ denotes the inner product of $f(\bm{x})$
and $g(\bm{x})$ in ${\cal L}^{2}({\cal X},\rho)$. By Mercer’s theorem
(Steinwart and Christmann, 2008), under regularity assumptions, the eigen-
expansion of $K(\bm{x},\bm{u})$ is $K(\bm{x},\bm{u})=\sum_{j\geq
1}\mu_{j}e_{j}(\bm{x})e_{j}(\bm{u}).$ Hence, the RKHS-norm of any $f\in{\cal
H}_{K}$ can also be expressed as

$\displaystyle\|f\|^{2}_{K}=\sum_{j\geq 1}\frac{\langle
f,e_{j}\rangle^{2}_{2}}{\mu_{j}}.$

The above result implies the decay rate of $\mu_{j}$ fully characterizes the
complexity of the RKHS.

The following technical assumptions are made to investigate the theoretical
properties of the AIPW in (3).

Assumption 1: There exists a positive constant $r\in(1/2,1]$, such that
$f^{*}(\bm{x})$ is in the range of the $r$th power of $L_{K}$, denoted as
$L_{K}^{r}$. Besides, the distribution of $\epsilon$ has a $q$-exponential
tail with some function $q(\cdot)$; that is, there exists a constant
$c_{1}>0$, such that $\Pr(\lvert\epsilon\rvert>t)\leq c_{1}\exp\\{-q(t)\\}$
for any $t>0$.

Assumption 2: There exist positive values $\kappa_{1,p}$ and $\kappa_{2,p}$,
which may depend on $p$, such that $\sup\limits_{\bm{x}\in{\cal
X}}\|K_{\bm{x}}\|_{K}\leq\kappa_{1,p}$ and $\sup\limits_{\bm{x}\in{\cal
X}}\|\partial_{l}K_{\bm{x}}\|_{K}\leq\kappa_{2,p}$ for $l=1,...,p$.

Assumption 3: There exists a positive constant $\xi_{1}>1$ such that

$\displaystyle\min_{l\in{\cal A}^{*}}\left\|g^{*}_{l}\right\|^{2}_{2}>$
$\displaystyle
c_{m}\max\left\\{\kappa_{1,p}\|f^{*}\|_{K},q^{-1}\left(\log\frac{4c_{1}m}{\delta_{m}}\right)\right\\}m^{-\frac{2r-1}{2(2r+1)}}(\log
p)^{\xi_{1}},$

where $c_{m}$ is provided in Lemma 1.

Assumption 4: There exists $\kappa\in(0,1/2)$ such that
$\kappa<\pi(\bm{x})<1-\kappa$ for all $\bm{x}\in{\cal X}$.

Assumption 5: $\mathbb{E}(\bm{X}\bm{X}^{\top})$ is invertible, and its
smallest eigenvalue is bounded away from zero by a fixed positive constant
$c_{min}$, and recall that $\bm{X}$ is the random vector associated with
$\bm{x}_{1},\ldots,\bm{x}_{n}$.

Assumption 6: Let $\bm{X}_{g}$ be the random vector associated with the $g$th
group, and we normalize $\bm{X}_{g}$ such that
$\mathbb{E}(\bm{X}_{g}^{\top}\bm{X}_{g})$ is a
$\mbox{df}_{g}\times\mbox{df}_{g}$ identity matrix. Then, there exists $L_{n}$
such that $\max_{\bm{x}}\max_{g}(\bm{x}_{g}^{\top}\bm{x}_{g})\leq nL_{n}^{2},$
where $\bm{x}_{g}$ corresponds to the normalized $\bm{X}_{g}$.

Assumption 7: $\max_{g=1,\ldots,G}\mbox{df}_{g}=O(1)$, there exists a constant
number $\zeta>0$ such that $\log(G)=o(n^{1/3-2\zeta})$ and $G\gg\log(n)$,
$N_{0}=O(1)$, $\lambda_{2}\asymp\log(G)$, i.e. $\lambda_{2}$ is of the order
$\log(G)$, and $L_{n}^{2}=O\\{1/\log(G)\\}$, where $N_{0}$ is the number of
non-zero group effects.

Assumptions 1–3 are proposed for the kernel-based sparse learning algorithm,
and Assumptions 4–7 are required by the group lasso logistic regression. In
Assumption 1, the integral operator $L_{K}$ is self-adjoint and semi-positive
definite, so its fractional operator $L_{K}^{r}$ is well-defined, and its
range is contained in ${\cal H}_{K}$ as long as $r\geq 1/2$; see Smale and
Zhou (2007) and Mendelson and Neeman (2010) for details. This implies that for
some function $h\in{\cal L}^{2}({\cal X},\rho)$, it holds
$L^{r}_{K}f^{*}=\sum_{j\geq 1}\mu^{r}_{j}\langle
h,e_{j}\rangle_{2}e_{j}\in{\cal H}_{K},$ so ensures strong estimation
consistency under the RKHS-norm. The second part of Assumption 1 characterizes
the tail behavior of the error distribution, and it relaxes the commonly-used
bounded response assumption in the machine learning literature (Smale and
Zhou, 2007; Rosasco et al., 2013; Lv et al., 2018). Besides, the assumption on
the error distribution is general and can be satisfied by a variety of
distributions (Wang and Leng, 2016; Zhang et al., 2016). For example, if the
error distribution is sub-Gaussian or bounded, then $q(t)=O(t^{2})$ suffices;
if the distribution of $\epsilon$ is sub-exponential,
$q(t)=O(\min\\{t/C,t^{2}/C^{2}\\})$ suffices with $C>0$. Assumption 2 bounds
the RKHS-norms associated with the kernel function and its gradient functions,
and it is satisfied by popular kernels, including the Gaussian kernel, linear
kernel and the Sobolev kernel (Smale and Zhou, 2007; Rosasco et al., 2013;
Yang et al., 2016). For example, $\kappa_{1,p}=\kappa_{2,p}=1$ for the
Gaussian kernel,
$K(\bm{x},\bm{u})=\exp\\{-\|\bm{x}-\bm{u}\|_{2}^{2}/(2\sigma^{2})\\}$, and
$\kappa_{1,p}=Cp$ and $\kappa_{2,p}=C$ for the linear kernel,
$K(\bm{x},\bm{u})=\bm{x}^{\top}\bm{u}$, for some positive constant $C$.
Assumption 3 requires that the gradient functions contain sufficient
information about the truly informative covariates. It is worthy pointing out
that we measure the significance of each gradient function to distinguish
informative and uninformative covariates without any explicit model
specification. The minimal signal strength in Assumption 3 is much tighter
than those in other nonparametric sparse learning methods (Huang et al., 2010;
Yang et al., 2016), which often require the signal strength to be bounded
below by some positive constant. Assumption 4 bounds the response probability,
and it is commonly used to avoid inefficient estimators. To obtain the desired
convergence rate, Assumption 4 guarantees that $m\asymp n$ in probability,
where $a_{n}\asymp b_{n}$ is equivalent to $a_{n}=O(b_{n})$ and
$b_{n}=O(a_{n})$. The smallest eigenvalue of $\mathbb{E}(\bm{X}\bm{X}^{\top})$
is bounded by a fixed positive constant in Assumption 5, and it is a special
case of assumption (b) of Meier et al. (2008). In Assumption 6, the
convergence rate of $L_{n}$ is related with that of the estimated response
probability. Assumption 7 is used to guarantee that

$\mathbb{E}\\{\lvert\eta_{\widehat{\bm{\beta}}_{\lambda_{2}}}(\bm{X})-\eta_{{\bm{\beta}}^{*}}(\bm{X})\rvert^{2}\\}=O_{\mathbb{P}}(n^{-2/3-2\zeta}),$
(9)

where $\eta_{\bm{\beta}}(\bm{x})=\beta_{0}+\bm{x}^{\top}\bm{\beta}_{1}$,
$\bm{\beta}^{*}=(\beta_{0}^{*},\bm{\beta}_{1}^{*})$ and the expectation is
taken with respect to $\bm{X}$ conditional on
$\widehat{\bm{\beta}}_{\lambda_{2}}$ or the observations. Specifically, the
value $\zeta$ is used to show $n^{-\zeta/2}\log(n)\to 0$ as $n\to 0$, so it
can be chosen arbitrarily small; see the proof of Theorem 2 for details.

######  Lemma 1.

Suppose Assumptions 1–2 are satisfied, and $\lambda={m^{-1/(2r+1)}}$. Then,
for any $\delta_{m}\in(0,1)$, with probability at least $1-\delta_{m}$, there
holds

$\displaystyle\max_{1\leq l\leq p}\
\big{\lvert}\|\widehat{g}_{l}\|^{2}_{m}-\|g^{*}_{l}\|_{2}^{2}\big{\rvert}\leq
c_{m}c_{p,q}{\log\Big{(}\frac{8p}{\delta_{m}}\Big{)}}{m^{-\Theta}},$

where $c_{m}$ is a constant depending only on $\kappa_{1,p},\kappa_{2,p}$ and
$\|f^{*}\|^{2}_{K}$,
$c_{p,q}=\max\left\\{\kappa_{1,p}\|f^{*}\|_{K},q^{-1}\left(\log\frac{4c_{1}m}{\delta_{m}}\right)\right\\}$
with $q^{-1}(\cdot)$ denoting the inverse function of $q(\cdot)$, and
$\Theta={\frac{2r-1}{2(2r+1)}}$.

Lemma 1 guarantees that $\left\|\widehat{g}_{l}\right\|^{2}_{m}$ converges to
$\|g_{l}^{*}\|_{2}^{2}$ with high probability, and it is crucial to establish
the selection consistency of the employed sparse learning algorithm. The
convergence result in Lemma 1 still holds even when the dimension diverges
with the sample size, and the quantities $\|f^{*}\|^{2}_{K}$ and
$\|L^{-r}_{K}f^{*}\|_{2}$, which may depend on the number of truly informative
covariates of $f^{*}(\bm{x})$, may also diverge as the sample size increases.
For instance, if $f^{*}(\bm{x})=\bm{x}^{T}\bm{\beta}^{*}$, then
$\|f^{*}\|^{2}_{K}=\|\bm{\beta}^{*}\|_{2}^{2}$, which clearly depends on the
number of truly informative covariates. However, such dependency is difficult
to quantify explicitly in a fully general case (Fukumizu and Leng, 2014).

The following lemma establishes the asymptotic selection consistency of the
proposed sparse learning method.

######  Lemma 2.

Suppose that the assumptions of Lemma 1 and Assumption 3 are satisfied. If
$v_{m}=0.5c_{m}c_{p,q}m^{-\Theta}(\log p)^{\xi_{1}},$ then
$\Pr\big{(}\widehat{\cal A}_{v_{m}}={\cal A}^{*}\big{)}\rightarrow 1,\ \
\mbox{as}\ \ m\rightarrow\infty.$

Lemma 2 shows that the selected covariates can exactly recover the truly
informative ones with probability tending to 1. This result is particularly
general in that it is established without any model specification. The
following theorem shows that $\widehat{f}_{0}(\bm{x})$ achieves a fast
convergence rate in term of the infinity norm, where $\widehat{f}_{0}(\bm{x})$
is obtained by the standard kernel ridge regression (6) based on the selected
covariates in $\widehat{\cal A}_{v_{n}}$.

######  Theorem 1.

Suppose the assumptions of Lemma 2 are satisfied and denote the probability
$\Pr(\widehat{\cal A}_{v_{n}}\neq{\cal A}^{*})=\Delta_{m}$. If
$\lambda={m^{-\frac{1}{2r+1}}}$, then with probability at least
$1-\delta_{m}-\Delta_{m}$, there holds

$\displaystyle\|\widehat{f}_{0}-f^{*}\|_{K}$ $\displaystyle\leq$
$\displaystyle
c_{m,2}c_{p_{0},q}{\log\Big{(}\frac{4}{\delta_{m}}\Big{)}}{m^{-\Theta}},$

where $p_{0}={|}{\cal A}^{*}{|}$,

$\displaystyle c_{m,2}$ $\displaystyle=$ $\displaystyle
4\max\\{\kappa^{2}_{2,p_{0}},\kappa^{2}_{2,p_{0}}\|f^{*}\|_{K},\|f^{*}\|^{2}_{K}\\}$
$\displaystyle\times\max\\{3\kappa_{1,p_{0}},2\sqrt{2}\kappa_{2,p_{0}}^{2},\|L_{K}^{-r}f^{*}\|_{2}\\},$

and
$c_{p_{0},q}=\max\left\\{\kappa_{1,p_{0}}\|f^{*}\|_{K},q^{-1}\left(\log\frac{4c_{1}m}{\delta_{m}}\right)\right\\}$.
Additionally, if we take $r=1$, and assume that $p_{0}=O(1)$, $\epsilon$ has
sub-Gaussian or sub-exponential tail and by Assumption 4, we have
$\left\|\widehat{f}_{0}-f^{*}\right\|_{\infty}=O_{\mathbb{P}}(n^{-\frac{1}{6}}\log
n).$

Theorem 1 establishes the uniform convergence rate of the refitted estimator
$\widehat{f}_{0}(\bm{x})$, and it plays a crucial role to establish the
central limit theory of the AIPW estimator in (3). The required tail behavior
of $\epsilon$ in Theorem 1 is to quantify $q^{-1}(\cdot)$ explicitly for
simplicity, and it can be extended to any error distribution satisfying
Assumption 1.

######  Lemma 3.

Suppose Assumptions 4–7 are satisfied. Then, given
$\widehat{\bm{\beta}}_{\lambda_{2}}$, there holds
$\mathbb{E}[\lvert\widehat{\pi}(\bm{X})-\pi^{*}(\bm{X})\rvert^{2}]=O_{\mathbb{P}}(n^{-2/3-2\zeta}).$

Lemma 3 establishes the convergence rate of the estimated response probability
using group lasso logistic regression. By Lemma 3, we essentially require that
the estimated response probability should be at least consistent. A similar
requirement is also discussed by Tan (2020). Specifically, Tan (2020) assumed
a correctly specified response model in order to achieve valid interval
estimator. If other penalized logistic regression estimators are considered,
Assumptions 4–7 should be replaced in order to guarantee Lemma 3.

By Theorem 1 and Lemma 3, we can validate the following central limit theorem
for the AIPW estimator in (3).

######  Theorem 2.

Suppose all the assumptions in Theorem 1 and Lemma 3 are satisfied. If
$\mathbb{E}\\{\lvert
f^{*}(\bm{X})+\delta\pi^{*}(\bm{X})^{-1}\\{Y-f^{*}(\bm{X})\\}\rvert^{2}\\}<\infty$,
then $\sqrt{n}(\widehat{\theta}_{AIPW}-\theta^{*})\to\rm{N}(0,\sigma^{2}),$ in
distribution, where $\pi^{*}(\bm{X})=\Pr(\delta=1\mid\bm{X})$ and
$\sigma^{2}=\mathop{\rm
var}{[f^{*}(\bm{X})+\delta\pi^{*}(\bm{X})^{-1}\\{Y-f^{*}(\bm{X})\\}]}$.

It is worthy pointing out that the derived result is particularly attractive
given the fact that the central limit theorem is built by nonparametric
estimation of $f^{*}(\bm{x})$ with diverging dimension, and to our knowledge,
such a result is novel in literature. More importantly, the variance term
$\sigma^{2}$ can be estimated by the sample variance of
$\\{\widehat{f}_{0}(\bm{x}_{i})+\delta_{i}\widehat{\pi}(\bm{x}_{i})^{-1}\\{y_{i}-\widehat{f}_{0}(\bm{x}_{i})\\}:i=1,\ldots,n\\}$:

$\widehat{\sigma}_{AIPW}^{2}=\frac{1}{n-1}\sum_{i=1}^{n}(\widehat{y}_{i}-\widehat{\theta}_{AIPW})^{2},$

where
$\widehat{y}_{i}=\widehat{f}_{0}(\bm{x}_{i})+\delta_{i}\widehat{\pi}(\bm{x}_{i})^{-1}\\{y_{i}-\widehat{f}_{0}(\bm{x}_{i})\\}$.
Thus, based on Theorem 2 and the estimated variance
$\widehat{\sigma}_{AIPW}^{2}$, we can also obtain the interval estimators of
$\theta^{*}$.

##  4 Numerical analysis

In this section, we compared the numerical performance of the proposed AIPW
estimator, denoted as Prop, against several state-of-the-art competitors under
two simulated experiments and a real-data application. For Prop, in all the
scenarios, we applied a Gaussian kernel,
$K(\bm{x},\bm{u})=\exp{\left(-{\|\bm{x}-\bm{u}\|_{2}^{2}}/(2\sigma_{n}^{2})\right)}$
with $\sigma_{n}$ being the median of all the pairwise distances among the
covariates (Jaakkola et al., 1999). As suggested by He et al. (2021), we also
applied the stability-based selection criterion (Sun et al., 2013) to
determine the thresholding value $v_{n}$ and set the ridge parameter
$\lambda_{n}=0.001$ for the employed sparse learning algorithm.

###  4.1 Simulated experiments

In this section, we considered $n\in\\{800,1\,000\\}$ and
$p\in\\{400,2\,000\\}$, and covariates were generated by
$x_{il}\sim\mbox{U}(-0.5,0.5)$ for $i=1,\ldots,n$ and $l=1,\ldots,p$, where
$x_{il}$ denoted the $l$th element of $\bm{x}_{i}$, and $\mbox{U}(-0.5,0.5)$
denoted a uniform distribution over $[-0.5,0.5]$. The following regression
models were applied to generate the response of interest:

  1. M1.

Linear regression model: $y_{i}=5x_{i1}+6x_{i2}+4x_{i3}+4x_{i4}+\epsilon_{i}$
with $\epsilon_{i}\sim\rm{N}(0,1)$.

  2. M2.

Nonlinear regression model:
$y_{i}=6x_{i1}+4(2x_{i2}+1)(2x_{i3}-1)+6h(x_{i4})+5\sin(x_{i5}\pi)/\\{2-\sin(x_{i5}\pi)\\}+\epsilon_{i}$,
where
$h(x)=0.1\sin(x_{i4}\pi)+0.2\cos(x_{i4}\pi)+0.3\sin(x_{i4}\pi)^{2}+0.4\cos(x_{i4}\pi)^{3}+0.5\sin(x_{i4}\pi)^{3}$,
and $\epsilon_{i}\sim\rm{N}(0,1)$.

For $i=1,\ldots,n$, the response indicator $\delta_{i}$ was generated by a
Bernoulli distribution with success probability $\pi^{*}(\bm{x}_{i})$, which
was obtained by the following models:

  1. R1.

Logistic response model:
$\mbox{logit}\\{\pi^{*}(\bm{x}_{i})\\}=-0.1+2x_{i1}+2x_{i3}$.

  2. R2.

Multi-modal response model $\pi^{*}(\bm{x}_{i})=\sin(6x_{i2}+8x_{i4})/3+0.5$.

The linear regression model M1 is commonly assumed in practice (Fan and Li,
2001). The nonlinear regression model M2, however, is more complex, and the
interaction effect is also taken into consideration. The logistic response
model R1 is widely used in practice. However, the response model R2 violates
(2.2), so it is used to test the robustness of the proposed AIPW estimator.

The primary interest was to estimate $\theta^{*}=\mathbb{E}(Y)$. For the
regression model M1, we had $\theta^{*}=0$. However, instead of deriving
$\theta^{*}$ analytically, we used $\tilde{\theta}=L^{-1}\sum_{l=1}^{L}y_{l}$
as the “true value” for the regression model M2, where
$\\{y_{l}:l=1,\ldots,L\\}$ was a random sample of size $L=1\,000\,000$. The
following competitors were considered:

  1. CC.

The sample mean of the complete cases,
$\widehat{\theta}_{cc}=m^{-1}\sum_{i=1}^{n}\delta_{i}y_{i}$, where
$m=\sum_{i=1}^{n}\delta_{i}$.

  2. PS.

Conventional propensity score estimator
$\widehat{\theta}_{ps}={n}^{-1}\sum_{i=1}^{n}\delta_{i}\pi^{-1}(\bm{x}_{i};\widehat{\bm{\beta}})y_{i},$
where
$\widehat{\bm{\beta}}^{\top}=(\widehat{\beta}_{0},\widehat{\bm{\beta}}_{1}^{\top})$
solves
$\sum_{i=1}^{n}\\{\delta_{i}-\pi(\bm{x}_{i};\widehat{\bm{\beta}})\\}(1,\bm{x}_{i}^{\top})=0$
without consideration of the sparsity.

  3. DI.

Deterministic imputation using kernel ridge regression (Wang and Kim, 2021)
$\widehat{\theta}_{di}={n}^{-1}\sum_{i=1}^{n}\\{\delta_{i}y_{i}+(1-\delta_{i})\widehat{f}(\bm{x}_{i})\\},$
where $\widehat{f}(\bm{x})$ is the fitted kernel ridge regression model based
on the fully observed data $\\{(\bm{x}_{i},y_{i}):\delta_{i}=1\\}$ without
employing sparse learning.

  4. NAIPW.

Naive AIPW estimator
$\widehat{\theta}_{AIPW1}={n}^{-1}\sum_{i=1}^{n}\\{\widehat{f}(\bm{x}_{i})+\delta_{i}\pi^{-1}(\bm{x}_{i};\widehat{\bm{\beta}})\\{y_{i}-\widehat{f}(\bm{x}_{i})\\}\\},$
where $\widehat{\bm{\beta}}$ is the same as that in the PS estimator without
consideration of the sparsity, and $\widehat{f}(\bm{x})$ is the same as that
in the DI estimator without employing sparse learning.

The CC estimator completely ignore the unobserved data, leading to a biased
estimator if $\mathbb{E}(\delta_{i}\mid\bm{x}_{i})$ involves covariates used
in the regression model. PS estimator is widely used in causal inference
(Rosenbaum and Rubin, 1983) and missing data analysis (Wooldridge, 2007). The
imputation methods are commonly used to provide a complete dataset, especially
in survey sampling; see Kim and Shao (2013, Chapter 4) for details. Recently,
Wang and Kim (2021) has proposed a kernel-based deterministic imputation
method, and we consider their method for comparison as well. Except for the
proposed AIPW estimator, we also considered the naive AIPW estimator based on
the conventional propensity score estimator and the deterministic imputation
estimator.

We conducted $M=500$ Monte Carlo simulations for each estimator under
different model setups. First, we compared different estimators in terms of
the Monte Carlo bias and the Monte Carlo standard error:

Bias $\displaystyle=\bar{\theta}_{n}^{(M)}-\theta,$ SE
$\displaystyle=\left\\{\frac{1}{M-1}\sum_{m=1}^{M}(\widehat{\theta}_{n}^{(m)}-\bar{\theta}_{n}^{(M)})^{2}\right\\}^{1/2},$

where $\bar{\theta}_{n}^{(M)}=M^{-1}\sum_{m=1}^{M}\widehat{\theta}_{n}^{(m)}$,
$\widehat{\theta}_{n}^{(m)}$ was a specific estimator of $\theta$ for the
$m$th Monte Carlo simulation. Simulation results were summarized in Table 1.
The CC estimator is biased since the response of interest $y_{i}$ is
correlated with the response index $\delta_{i}$. Even though a logistic model
is correctly specified for the response model R1, the PS estimator is still
biased or even unrealistic due to the curse of dimensionality. Since the NAIPW
estimator was obtained using the same response model as the PS estimator, it
is also questionable, especially when the sample size is small and the number
of useless covariates is large. Since the response probability was not used by
the DI estimator, it does not suffer the same problem as the PS estimator.
However, even under the linear regression model, the bias of the DI estimator
may not be negligible compared with its standard error. Compared with its
competitors, the proposed AIPW estimator performs the best since it has the
smallest bias under most model setups, and its standard error is reasonably
small.

Table 1: Summary of the Monte Carlo bias (Bias) and standard error (SE) corresponding to the five estimators under different model setups, and the unit is 0.1. For “Model”, “C1–C4” represent (M1, R1),(M2, R1), (M2, R2) and (M2, R2), respectively. For “Size”, “I–IV” corresponds to $(n,p)=(800,400)$, $(1\,000,400)$, $(800,2\,000)$, and $(1\,000,2\,000)$, respectively. Notation “-” is used when the absolute value of either bias or standard error is greater than $100$.  Model | Size | CC | PS | DI | NAIPW | Prop
---|---|---|---|---|---|---
| I | 6.9 (1.5) | \- (-) | 0.9 (1.1) | \- (-) | 0.3 (1.2)
| II | 7.0 (1.3) | \- (-) | 0.9 (1.0) | \- (-) | 0.1 (1.0)
| III | 7.0 (1.5) | 3.4 (0.7) | 1.1 (1.2) | 1.1 (1.2) | 0.4 (1.3)
C1 | IV | 6.9 (1.3) | 3.3 (0.6) | 1.0 (1.0) | 1.0 (1.0) | 0.2 (1.1)
| I | -0.7 (2.8) | \- (-) | 0.1 (2.0) | \- (-) | 0.0 (2.0)
| II | -0.8 (2.4) | \- (-) | -0.1 (1.7) | \- (-) | -0.1 (1.6)
| III | -0.7 (2.6) | 2.5 (1.3) | 0.3 (2.1) | 0.3 (2.1) | 0.0 (2.0)
C2 | IV | -0.6 (2.3) | 2.6 (1.1) | 0.0 (1.7) | 0.0 (1.7) | -0.1 (1.7)
| I | -1.2 (1.5) | \- (-) | -0.3 (1.1) | \- (-) | -0.2 (1.1)
| II | -1.3 (1.3) | 0.3 (57.1) | -0.3 (1.0) | 0.2 (8.4) | -0.2 (1.0)
| III | -1.2 (1.5) | -0.6 (0.7) | -0.2 (1.1) | -0.2 (1.1) | -0.1 (1.1)
C3 | IV | -1.2 (1.3) | -0.6 (0.7) | -0.2 (1.0) | -0.2 (1.0) | -0.2 (1.0)
| I | 1.7 (2.6) | \- (-) | 0.2 (1.9) | \- (-) | 0.0 (1.9)
| II | 1.7 (83.6) | -3.0 (2.2) | 0.2 (1.6) | -0.2 (17.9) | 0.0 (1.7)
| III | 1.8 (2.6) | 3.6 (1.3) | 0.3 (2.0) | 0.3 (2.0) | 0.1 (2.0)
C4 | IV | 1.9 (2.2) | 3.6 (1.1) | 0.3 (1.6) | 0.3 (1.6) | 0.0 (1.6)

Next, the proposed AIPW estimator was evaluated by the relative bias of the
variance estimator and its coverage rate of a 95% confidence interval:

RB
$\displaystyle=\frac{\bar{\sigma}_{n}^{2(M)}-\mbox{SE}^{2}}{\mbox{SE}^{2}},$
CR
$\displaystyle=\frac{1}{M}\sum_{m=1}^{M}\mathbb{I}(\widehat{\theta}_{n}^{(m)}-1.96\widehat{\sigma}_{n}^{(m)}\leq\theta\leq\widehat{\theta}_{n}^{(m)}+1.96\widehat{\sigma}_{n}^{(m)}),$

where
$\bar{\sigma}_{n}^{2(M)}=M^{-1}\sum_{i=1}^{M}\widehat{\sigma}_{n}^{2(m)}$,
$\widehat{\sigma}_{n}^{2(m)}$ is the variance estimator for the $m$th Monte
Carlo simulation, $\widehat{\sigma}_{n}^{(m)}$ is the square root of
$\widehat{\sigma}_{n}^{2(m)}$, and $\mathbb{I}(a\leq x\leq b)$ is an indicator
function of $u$ for given $a\leq b$, and it takes value 1 if $u\in[a,b]$ and 0
otherwise. We conducted 500 Monte Carlo simulations, and Table 2 summarized
the corresponding results. Generally, as the sample size increases from
$n=800$ to $n=1\,000$, the relative bias of the variance estimator decreases,
and it is negligible if sample size $n=1\,000$. Thus, the proposed variance
estimator performs well, especially when the sample size is large. The
coverage rates are close to their nominal truth 0.95 when sample size is
large. Since the variance of the proposed AIPW estimator is under-estimated
for the setup with regression model $M1$ and response model $R1$, the
corresponding coverage rate is much lower than 0.95. As the sample size
increases, however, the coverage rate gets closer to its nominal truth. For
the two setups with response model R2, a logistic regression model is wrongly
specified for the response indicator. However, the absolute values of the
relative bias of the variance estimator are generally less than 0.05 and the
corresponding coverage rates are close to the nominal truth 0.95, specially
when the sample size is large. Thus, the proposed AIPW estimator is indeed
robust against a wrongly specified response model.

Table 2: Relative bias (RB) of the variance estimator and coverage rate (CR) of the 95% confidence interval of $\theta$ for the proposed method. For “Model”, “C1–C4” represent (M1, R1),(M2, R1), (M2, R2) and (M2, R2), respectively. “I–IV” corresponds to $(n,p)=(800,400)$, $(1\,000,400)$, $(800,2\,000)$, and $(1\,000,2\,000)$, respectively.  | RB |  | CR
---|---|---|---
Model | I | II | III | IV |  | I | II | III | IV
C1 | -0.10 | -0.07 | -0.33 | -0.20 |  | 0.92 | 0.94 | 0.88 | 0.92
C2 | -0.13 | 0.01 | -0.14 | -0.01 |  | 0.94 | 0.94 | 0.93 | 0.95
C3 | 0.04 | -0.02 | -0.03 | 0.05 |  | 0.95 | 0.95 | 0.95 | 0.96
C4 | -0.10 | -0.01 | -0.14 | 0.05 |  | 0.94 | 0.94 | 0.93 | 0.96

###  4.2 Application to a supermarket dataset

In this section, the proposed AIPW estimator and its competitors were applied
to a supermarket dataset (Wang, 2009), which was collected from a major
supermarket located in northern China, consisting of daily sale records of
$p=6\,398$ products on $n=464$ days. This data included almost all kinds of
daily necessities and the response of interest was the number of customers on
each day, and the covariates are the daily sale volumes of each product. For
simplicity, denote $y_{i}$ and $\bm{x}_{i}=(x_{i1},\ldots,x_{ip})^{\top}$ to
the be the response of interest and the corresponding covariate for the $i$th
day. In this section, we were interested in estimating the average number of
customers visiting the supermarket. This dataset was fully observed, and it
was studentized before analyzing. Thus, the sample mean
$\widehat{\theta}=n^{-1}\sum_{i=1}^{n}y_{i}=0$ served as a benchmark.

To compare the performance of the proposed AIPW estimator with other
competitors, we considered the following missing mechanism for $y_{i}$:

$\displaystyle\mbox{logit}\\{\pi(\bm{x}_{i})\\}=1-0.6x_{i5}-x_{i6}+0.5x_{i10},$
(10)

and $y_{i}$ was treated as observed if and only if $\delta_{i}=1$, where
$\delta_{i}=1$ with probability $\pi_{i}$. The mechanism in (10) was MAR, and
the corresponding covariates were identified as informative for estimating the
response of interest by He et al. (2021). Then, instead of observing the whole
data, we assumed that all the covariates and only $\\{y_{i}:\delta_{i}=1\\}$
were available, and the resulting response rate was about 0.70.

We generated 500 incomplete datasets using (10) and compare the estimators in
Section 4.1. Table 3 summarized the average of the estimators and the
corresponding standard error. The CC estimator is highly biased since it
ignores the missing mechanism. The performance of the PS estimator is also
questionable in that the response model using all covariates results in
overfitting. The DI, NAIPW and Prop estimators outperform the CC and PS
estimators since their estimates are much closer to 0. However, the standard
error of the Prop estimator is much smaller than the other two, illustrating
the superior of the proposed AIPW estimator.

Table 3: The average and standard error of 500 incomplete datasets for estimating the number of customers visiting the supermarket. Since studentized is applied, the sample mean 0 serves as the benchmark.  | CC | PS | DI | NAIPW | Prop
---|---|---|---|---|---
Estimate | -0.20 | -0.13 | -0.04 | -0.04 | -0.04
Standard error | 0.03 | 0.02 | 0.03 | 0.03 | 0.01

##  5 Conclusion

In this paper, we propose a novel AIPW estimator to infer the population mean,
which incorporates an efficient nonparametric imputation with sparse structure
and a penalized propensity score estimator under the assumption of missing at
random. The proposed method is computationally efficient and allows the
dimension diverging. More importantly, the estimation consistency as well as
the corresponding central limit theorem are established under regularity
assumptions. Its superior is also supported by several simulated examples and
one application to a supermarket dataset.

## Acknowledgement

Xin He’s research is supported by NSFC Grant No. 11901375; Xiaojun Mao’s
research is supported by NSFC Grant No. 12001109 and 92046021, the Science and
Technology Commission of Shanghai Municipality grant 20dz1200600. Zhonglei
Wang’s research is supported by NSFC Grant No. 11901487 and 72033002,
Fundamental Scientific Center of National Natural Science Foundation of China
Grant No. 71988101.

##  Appendix A Block co-ordinate gradient descent algorithm

For $g$th group of $\bm{\beta}_{1}$, consider a vector $\bm{d}$ such that
$\bm{d}_{k}=0$ for $k\neq g$, and assume that the
$\mbox{df}_{g}\times\mbox{df}_{g}$ submatrix is of the form
$H_{gg}^{(t)}=h_{g}^{(t)}I_{\scriptsize\mbox{df}_{g}}$ for some scalar
$h_{g}^{(t)}$, where $I_{m}$ is an $m\times m$ identity matrix.

If $\lVert\nabla
l(\widehat{\bm{\beta}}^{(t)})_{g}-h_{g}^{(t)}\widehat{\bm{\beta}}^{(t)}_{g}\rVert_{2}\leq\lambda_{2}\mbox{df}_{g}^{1/2}$,
let $\bm{d}_{g}^{(t)}=-\widehat{\bm{\beta}}^{(t)}_{g}$. Otherwise,

$\bm{d}_{g}^{(t)}=-\frac{1}{h_{g}^{(t)}}\left\\{\nabla
l(\widehat{\bm{\beta}}^{(t)})_{g}-\lambda_{2}\mbox{df}_{g}^{1/2}\frac{\nabla
l(\widehat{\bm{\beta}}^{(t)})_{g}-h_{g}^{(t)}\widehat{\bm{\beta}}^{(t)}_{g}}{\lVert\nabla
l(\widehat{\bm{\beta}}^{(t)})_{g}-h_{g}^{(t)}\widehat{\bm{\beta}}^{(t)}_{g}\rVert_{2}}\right\\}.$

If $\bm{d}^{(t)}\neq 0$, let
$\bm{\beta}^{(t+1)}=\bm{\beta}^{(t)}+\alpha^{(t)}\bm{d}^{(t)}$, where
$\alpha^{(t)}$ is the largest value among $\\{\alpha_{0}\delta^{l}:l\geq 0\\}$
such that

$S_{\lambda_{2}}(\bm{\beta}^{(t)}+\alpha^{(t)}\bm{d}^{(t)})-S_{\lambda_{2}}(\bm{\beta}^{(t)})\leq\alpha^{(t)}\sigma\Delta^{(t)},$

$\delta\in(0,1)$, $\sigma\in(0,1)$, $\alpha_{0}>0$, and

$\Delta^{(t)}=-\left(\bm{d}^{(t)}\right)^{\top}\nabla
l(\widehat{\bm{\beta}}^{(t)})+\lambda_{2}\mbox{df}_{g}^{1/2}\lVert\widehat{\bm{\beta}}^{(t)}_{g}+\bm{d}^{(t)}_{g}\rVert_{2}-\lambda_{2}\mbox{df}_{g}^{1/2}\lVert\widehat{\bm{\beta}}^{(t)}_{g}\rVert_{2}.$

See Meier et al. (2008) for details.

##  Appendix B Proofs

######  Proposition 1.

Suppose Assumptions 1–2 are satisfied. Then, with probability at least
$1-\delta_{n}/2$, there holds

$\displaystyle\big{\|}\widehat{f}-f^{*}\big{\|}_{K}\leq
2\log\left(\frac{8}{\delta_{n}}\right)\Big{[}\frac{3\kappa_{1,p}}{n^{1/2}\lambda_{n}}\left\\{\kappa_{1,p}\|f^{*}\|_{K}+q^{-1}\left(\log\frac{4c_{1}n}{\delta_{n}}\right)\right\\}$
$\displaystyle+\lambda_{n}^{r-{1}/{2}}\|L^{-r}_{K}f^{*}\|_{2}\Big{]}.$

The proof of Proposition 1 is similar as that in He et al. (2021) and thus we
omit it here.

###### Proof of Lemma 1.

The proof of Lemma 1 is similar as that in He et al. (2021) by using
Proposition 1, the property of Hilbert-Schmidt operators and the concentration
inequalities in Hilbert-Schmidt operator space. Thus we omit the detail here.
∎

###### Proof of Lemma 2.

The proof of Lemma 2 is similar as that in He et al. (2021), and thus we omit
the detail here. ∎

###### Proof of Theorem 1.

Define the event that

$\displaystyle{\cal
C}_{1}=\left\\{\|\widehat{f}_{0}-f^{*}\|_{\infty}>c_{m,2}\max\left\\{\kappa_{1,p_{0}}\|f^{*}\|_{K},q^{-1}\left(\log\frac{2c_{1}m}{\delta_{m}}\right)\right\\}\right.$
$\displaystyle\times\left.{\log\left(\frac{4}{\delta_{m}}\right)}{m^{-\frac{2r-1}{2(2r+1)}}}\right\\}.$
(11)

Then, the probability $\Pr({\cal C}_{1})$ can be decomposed as

$\displaystyle\Pr\left({\cal C}_{1}\right)$ $\displaystyle=\Pr\left({\cal
C}_{1},\big{\\{}\widehat{\cal A}_{v_{m}}={\cal
A}^{*}\big{\\}}\right)+\Pr\left({\cal C}_{1},\big{\\{}\widehat{\cal
A}_{v_{m}}={\cal A}^{*}\big{\\}}\right)$ $\displaystyle=\Pr\left({\cal
C}_{1}\mid\big{\\{}\widehat{\cal A}_{v_{m}}={\cal
A}^{*}\big{\\}}\right)\Pr\left(\widehat{\cal A}_{v_{m}}={\cal A}^{*}\right)$
$\displaystyle\quad+\Pr\left({\cal C}_{1}\mid\big{\\{}\widehat{\cal
A}_{v_{m}}={\cal A}^{*}\big{\\}}\right)\Pr\left(\widehat{\cal
A}_{v_{m}}\neq{\cal A}^{*}\right)$ $\displaystyle\leq\Pr\left({\cal
C}_{1}\mid\big{\\{}\widehat{\cal A}_{v_{m}}={\cal
A}^{*}\big{\\}}\right)(1-\Delta_{m})+\Delta_{m}.$

By Lemma 2, we have $\Delta_{m}\rightarrow 0$ and $(1-\Delta_{m})\rightarrow
1$. For $\Pr({\cal C}_{1}\mid\big{\\{}\widehat{\cal A}_{v_{m}}={\cal
A}^{*}\big{\\}})$, by applying the proof in Proposition 1 conditioning on
$\big{\\{}\widehat{\cal A}_{v_{m}}={\cal A}^{*}\big{\\}}$, with probability at
least $1-\delta_{m}$, there holds

$\displaystyle\big{\|}\widehat{f}-f^{*}\big{\|}_{K}$ $\displaystyle\leq
c_{m,2}\max\left\\{\kappa_{1,p_{0}}\|f^{*}\|_{K},q^{-1}(\log\frac{2c_{1}m}{\delta_{m}})\right\\}$
$\displaystyle\times{\log\left(\frac{4}{\delta_{m}}\right)}{m^{-\frac{2r-1}{2(2r+1)}}},$

which implies $\Pr({\cal C}_{1}\mid\big{\\{}\widehat{\cal A}_{v_{m}}={\cal
A}^{*}\big{\\}})\leq\delta_{m}$. Combining the above results, we have
$\Pr({\cal C}_{1})\leq\delta_{m}+\Delta_{m}$. This completes the proof of the
first part in Theorem 1.

Additionally, by Assumption 4, we have $m=O(n)$, and if we take $r=1$ and
assume that $p_{0}=O(1)$ and $\epsilon$ has sub-Gaussian or sub-exponential
tail, there holds

$\big{\|}\widehat{f}-f^{*}\big{\|}_{K}=O_{\mathbb{P}}(n^{-\frac{1}{6}}\log
n).$

Note that

$\displaystyle\big{\|}\widehat{f}-f^{*}\big{\|}_{\infty}$ $\displaystyle=$
$\displaystyle\sup_{\bm{x}}\lvert\widehat{f}(\bm{x})-f^{*}(\bm{x})\rvert$
$\displaystyle=$
$\displaystyle\sup_{\bm{x}}\lvert\langle\widehat{f}-f^{*},K_{\bm{x}}\rangle_{K}\rvert\leq\kappa_{1,p_{0}}\big{\|}\widehat{f}-f^{*}\big{\|}_{K},$

which completes the proof. ∎

###### Proof of Lemma 3.

By Assumptions 4–7, Meier et al. (2008) showed (9). Denote
$g(x)=\\{1+\exp(-x)\\}^{-1}$, and we can show that
$g^{\prime}(x)={\mbox{d}g(x)}/{\mbox{d}x}=g(x)\\{1-g(x)\\}.$ That is,

$\lvert g^{\prime}(x)\rvert\leq 1,$ (12)

for any $x$ by the fact that $0\leq g(x)\leq 1$. Thus, by (12) and the mean
value theorem, we conclude that $g(x)$ is Lipschitz continuous in the sense
that

$\lvert g(x_{1})-g(x_{2})\rvert\leq\lvert x_{1}-x_{2}\rvert,$

for any $x_{1}$ and $x_{2}$ in $\mathbb{R}$. By noting the fact that
$\widehat{\pi}(\bm{x})=g\\{\eta_{\widehat{\bm{\beta}}_{\lambda_{2}}}(\bm{x})\\}$
and $\pi^{*}(\bm{x})=g\\{\eta_{\bm{\beta}_{0}}(\bm{x})\\}$, by (9), we have

$\displaystyle\mathbb{E}\\{\lvert\widehat{\pi}(\bm{X})-\pi^{*}(\bm{X})\rvert^{2}\\}$
$\displaystyle\leq$
$\displaystyle\mathbb{E}\\{\lvert\eta_{\widehat{\bm{\beta}}_{\lambda_{2}}}(\bm{X})-\eta_{{\bm{\beta}}_{0}}(\bm{X})\rvert^{2}\\}$
(13) $\displaystyle=$ $\displaystyle O_{\mathbb{P}}(n^{-2/3-2\zeta}),$

where the expectation is taken conditional on
$\widehat{\bm{\beta}}_{\lambda_{2}}$. By (13), we have shown Lemma 3.

∎

######  Lemma 4.

Suppose Assumptions 4–7 are satisfied. Then, given
$\widehat{\bm{\beta}}_{\lambda_{2}}$, there holds

$\max\\{\lvert\widehat{\beta}_{k}-\beta_{k}^{*}\rvert:k=0,\ldots,p\\}=O_{\mathbb{P}}(n^{-1/3-\zeta}),$

where $\widehat{\beta}_{k}$ and $\beta_{k}^{*}$ are the $(k+1)$th component of
$\widehat{\bm{\beta}}_{\lambda_{2}}$ and $\bm{\beta}^{*}$, respectively.

###### Proof of Lemma 4.

Given the estimated parameters for the response model, (9) can be re-expressed
as

$\displaystyle\mathbb{E}\\{\lvert\eta_{\widehat{\bm{\beta}}_{\lambda_{2}}}(\bm{X})-\eta_{{\bm{\beta}}^{*}}(\bm{X})\rvert^{2}\\}$
(14) $\displaystyle=$
$\displaystyle\mathbb{E}\\{(\widehat{\bm{\beta}}_{\lambda_{2}}-\bm{\beta}^{*})^{\top}\bm{X}\bm{X}^{\top}(\widehat{\bm{\beta}}_{\lambda_{2}}-\bm{\beta}^{*})\\}$
$\displaystyle=$
$\displaystyle(\widehat{\bm{\beta}}_{\lambda_{2}}-\bm{\beta}^{*})^{\top}\mathbb{E}(\bm{X}\bm{X}^{\top})(\widehat{\bm{\beta}}_{\lambda_{2}}-\bm{\beta}^{*})\\}$
$\displaystyle=$ $\displaystyle O_{\mathbb{P}}(n^{-2/3-2\zeta}),$

By Assumption A5 and (14), we conclude that

$\displaystyle(\widehat{\bm{\beta}}_{\lambda_{2}}-\bm{\beta}^{*})^{\top}(\widehat{\bm{\beta}}_{\lambda_{2}}-\bm{\beta}^{*})$
$\displaystyle=$
$\displaystyle(\widehat{\beta}_{0}-\beta_{0}^{*})^{2}+\sum_{k=1}^{p}(\widehat{\beta}_{k}-\beta_{k}^{*})^{2}$
(15) $\displaystyle=$ $\displaystyle O_{\mathbb{P}}(n^{-2/3-2\zeta}).$

Notice that

$\max\\{(\widehat{\beta}_{k}-\beta_{k}^{*})^{2}:k=0,\ldots,p\\}\leq(\widehat{\beta}_{0}-\beta_{0}^{*})^{2}+\sum_{k=1}^{p}(\widehat{\beta}_{k}-\beta_{k}^{*})^{2}.$
(16)

Thus, we have proved Lemma 4 by (15) and (16). ∎

###### Proof of Theorem 2.

For simplicity, denote $\pi^{*}_{i}=\pi^{*}(\bm{x}_{i})$ and
$\widehat{\pi}_{i}=\widehat{\pi}(\bm{x}_{i})=(1+\exp[-\\{(1,\bm{x}_{i}^{\top})\widehat{\bm{\beta}}_{\lambda_{2}}\\}])^{-1}$.
What if we consider

$\widehat{\theta}_{AIPW}=\frac{1}{n}\sum_{i=1}^{n}\left[\widehat{f}_{0}(\bm{x}_{i})+\frac{\delta_{i}}{\widehat{\pi}_{i}}\left\\{y_{i}-\widehat{f}_{0}(\bm{x}_{i})\right\\}\right],$

where $\widehat{\pi}_{i}$ is an estimator of $\Pr(\delta_{i}=1\mid\bm{x}_{i})$
by the group lasso for logistic regression.

Then, we have

$\displaystyle\widehat{\theta}_{AIPW}$ $\displaystyle=$
$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\left[{f}^{*}(\bm{x}_{i})+\left\\{\widehat{f}_{0}(\bm{x}_{i})-{f}^{*}(\bm{x}_{i})\right\\}+\frac{\delta_{i}}{{\pi}_{i}^{*}}\\{y_{i}-{f}^{*}(\bm{x}_{i})\\}\right.$
(17)
$\displaystyle+\frac{\delta_{i}}{\widehat{\pi}_{i}}\\{y_{i}-{f}^{*}(\bm{x}_{i})\\}-\frac{\delta_{i}}{{\pi}_{i}^{*}}\\{y_{i}-{f}^{*}(\bm{x}_{i})\\}$
$\displaystyle\left.+\frac{\delta_{i}}{\widehat{\pi}_{i}}\\{{f}^{*}(\bm{x}_{i})-\widehat{f}_{0}(\bm{x}_{i})\\}\right]$
$\displaystyle=$
$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\left[{f}^{*}(\bm{x}_{i})+\frac{\delta_{i}}{{\pi}_{i}^{*}}\\{y_{i}-{f}^{*}(\bm{x}_{i})\\}\right]$
$\displaystyle+\frac{1}{n}\sum_{i=1}^{n}\left\\{1-\frac{\delta_{i}}{\widehat{\pi}_{i}}\right\\}\left\\{\widehat{f}_{0}(\bm{x}_{i})-{f}^{*}(\bm{x}_{i})\right\\}$
$\displaystyle+\frac{1}{n}\sum_{i=1}^{n}\left[\frac{\delta_{i}}{\widehat{\pi}_{i}}-\frac{\delta_{i}}{{\pi}_{i}^{*}}\right]\\{y_{i}-{f}^{*}(\bm{x}_{i})\\}$
$\displaystyle=$
$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\left[{f}^{*}(\bm{x}_{i})+\frac{\delta_{i}}{{\pi}_{i}^{*}}\\{y_{i}-{f}^{*}(\bm{x}_{i})\\}\right]$
$\displaystyle+\frac{1}{n}\sum_{i=1}^{n}\left\\{1-\frac{\delta_{i}}{\widehat{\pi}_{i}}\right\\}\left\\{\widehat{f}_{0}(\bm{x}_{i})-{f}^{*}(\bm{x}_{i})\right\\}$
$\displaystyle+\frac{1}{n}\sum_{i=1}^{n}\left[\frac{\delta_{i}}{\widehat{\pi}_{i}}-\frac{\delta_{i}}{{\pi}_{i}^{*}}\right]\epsilon_{i},$

where $\epsilon_{i}=y_{i}-{f}^{*}(\bm{x}_{i})$.

First, we consider the first term of (17), and we have

$\displaystyle\mathbb{E}\left[f^{*}(\bm{X})+\frac{\delta}{\pi^{*}(\bm{X})}\\{Y-f^{*}(\bm{X})\\}\right]$
$\displaystyle=$
$\displaystyle\mathbb{E}\left(\mathbb{E}\left[f^{*}(\bm{X})+\frac{\delta}{\pi^{*}(\bm{X})}\\{Y-f^{*}(\bm{X})\\}\right]\mid\bm{X},Y\right)$
$\displaystyle=$ $\displaystyle\mathbb{E}(Y),$

where $\delta$ is a binary random variable with success probability
$\pi^{*}(\bm{X})$ conditional on $\bm{X}$. Since $\mathbb{E}\\{\lvert
f^{*}(\bm{X})+\delta\pi^{*}(\bm{X})^{-1}\\{Y-f^{*}(\bm{X})\\}\rvert^{2}\\}<\infty$,
by the classical central limit theorem (Van der Vaart, 2000, Example 2.1), we
have

$\sqrt{n}\left(\frac{1}{n}\sum_{i=1}^{n}\left[{f}^{*}(\bm{x}_{i})+\frac{\delta_{i}}{{\pi}_{i}^{*}}\\{y_{i}-{f}^{*}(\bm{x}_{i})\\}\right]-\theta^{*}\right)\to\rm{N}(0,\sigma^{2}),$
(18)

in distribution under regularity conditions, where $\theta=\mathbb{E}(y)$ and
$\sigma^{2}$ is to be estimated.

Next, we consider the third term of (17) .

$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\left[\frac{\delta_{i}}{\widehat{\pi}_{i}}-\frac{\delta_{i}}{{\pi}_{i}^{*}}\right]\epsilon_{i}$
$\displaystyle=$
$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\frac{\delta_{i}\epsilon_{i}}{{\pi}_{i}^{*}}\frac{{\pi}_{i}^{*}-\widehat{\pi}_{i}}{\widehat{\pi}_{i}}.$

By Assumption 4 and Lemma 3, we conclude that
$({\pi}_{i}^{*}-\widehat{\pi}_{i})\widehat{\pi}_{i}^{-1}=o_{\mathbb{P}}(1)$
uniformly for $i=1,\ldots,n$. Since

$\frac{1}{n}\sum_{i=1}^{n}\frac{\delta_{i}\epsilon_{i}}{{\pi}_{i}^{*}}=O_{\mathbb{P}}(n^{-1/2}),$

we conclude that

$\frac{1}{n}\sum_{i=1}^{n}\left[\frac{\delta_{i}}{\widehat{\pi}_{i}}-\frac{\delta_{i}}{{\pi}_{i}^{*}}\right]\epsilon_{i}=o_{\mathbb{P}}(n^{-1/2}).$
(19)

Now, we focus on the second term of (17), and consider

$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\left\\{1-\frac{\delta_{i}}{\widehat{\pi}_{i}}\right\\}\left\\{\widehat{f}_{0}(\bm{x}_{i})-{f}^{*}(\bm{x}_{i})\right\\}$
$\displaystyle=$
$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\left\\{1-\frac{\delta_{i}}{{\pi}_{i}^{*}}\right\\}\left\\{\widehat{f}_{0}(\bm{x}_{i})-{f}^{*}(\bm{x}_{i})\right\\}$
$\displaystyle+\frac{1}{n}\sum_{i=1}^{n}\frac{\delta_{i}}{{{\pi}_{i}^{*}}\widehat{\pi}_{i}}(\widehat{\pi}_{i}-{\pi}_{i}^{*})\left\\{\widehat{f}_{0}(\bm{x}_{i})-{f}^{*}(\bm{x}_{i})\right\\}.$

By Theorem 1,
$\widehat{f}_{0}(\bm{x}_{i})-{f}^{*}(\bm{x}_{i})=O_{\mathbb{P}}(\log(n)n^{-1/6})$
uniformly for $i=1,\ldots,n$. Since
$n^{-1}\sum_{i=1}^{n}\\{1-\delta_{i}(\pi_{i}^{*})^{-1}\\}=O_{\mathbb{P}}(n^{-1/2})$,
the first term of is of (LABEL:eq:_A3) is of the order
$o_{\mathbb{P}}(n^{-1/2})$. Besides, to show the second term of (LABEL:eq:_A3)
is also of the order $o_{\mathbb{P}}(n^{-1/2})$, it is enough to show

$\frac{1}{n}\sum_{i=1}^{n}\frac{\delta_{i}}{{{\pi}_{i}^{*}}\widehat{\pi}_{i}}(\widehat{\pi}_{i}-{\pi}_{i}^{*})=O_{\mathbb{P}}(n^{-1/3-\zeta/2}),$
(21)

where $\zeta$ is in Assumption 7.

By Lemma 4, we conclude that
$\max\\{\lvert\widehat{\beta}_{k}-\beta_{k}^{*}\rvert:k=0,\ldots,p\\}=o_{p}(1).$
Denote $A_{n}$ to be the event that
$\\{\max\\{\lvert\widehat{\beta}_{k}-\beta_{k}^{*}\rvert:k=0,\ldots,p\\}\geq
C_{\kappa}\\}$, where $C_{\kappa}$ is a positive constant such that
$\min\\{\widehat{\bm{\beta}}_{\lambda_{2}}(\bm{x})\geq\kappa/2:\bm{x}\in\mathcal{X}\\}$.
The existence of $C_{\kappa}$ is guaranteed by the compactness of
$\mathcal{X}$ and Assumption A4. Then, we have $P(A_{n})\to 0$ as
$n\to\infty$. On $\bm{x}\in A_{n}^{C}$, we conclude
$\hat{\pi}(\bm{x})\geq\kappa/2$.

Since $\\{\bm{x}_{1},\ldots,\bm{x}_{n}\\}$ is a random sample, given
$\widehat{\bm{\beta}}_{\lambda_{2}}$, for any positive constant $C$, we
consider

$\displaystyle\Pr\left(\left\lvert\frac{1}{n}\sum_{i=1}^{n}\frac{\delta_{i}}{\pi^{*}(\bm{X}_{i})\widehat{\pi}(\bm{X}_{i})}\\{\widehat{\pi}(\bm{X}_{i})-\pi^{*}(\bm{X}_{i})\\}\right\rvert\geq
Cn^{-1/3-\zeta/2}\right)$ (22) $\displaystyle\leq$ $\displaystyle
Pr\left(\left\lvert\frac{1}{n}\sum_{i=1}^{n}\frac{2\delta_{i}}{\kappa^{2}}\\{\widehat{\pi}(\bm{X}_{i})-\pi^{*}(\bm{X}_{i})\\}\right\rvert\geq
Cn^{-1/3-\zeta/2}\right)+P(A_{n})$ $\displaystyle\leq$
$\displaystyle\frac{\mathbb{E}[2\kappa^{-2}n^{-1}\sum_{i=1}^{n}[\delta_{i}\\{\pi^{*}(\bm{X}_{i})\widehat{\pi}(\bm{X}_{i})\\}^{-1}\\{\widehat{\pi}(\bm{X}_{i})-\pi^{*}(\bm{X}_{i})\\}]^{2}}{C^{2}n^{-2/3-\zeta}}$
$\displaystyle+P(A_{n})$ $\displaystyle\leq$
$\displaystyle\frac{2\sum_{i=1}^{n}\mathbb{E}\\{\widehat{\pi}(\bm{X}_{i})-\pi^{*}(\bm{X}_{i})\\}^{2}}{n\kappa^{2}C^{2}n^{-2/3-\zeta}}+P(A_{n})$
$\displaystyle\leq$
$\displaystyle\frac{\mathbb{E}\\{\widehat{\pi}(\bm{X}_{1})-\pi^{*}(\bm{X}_{1})\\}^{2}}{\kappa^{2}C^{2}n^{-2/3-\zeta}}+P(A_{n})$
$\displaystyle=$ $\displaystyle o_{\mathbb{P}}(1),$

where $\bm{X}_{i}$ is the random variable associated with $\bm{x}_{i}$, the
first inequality holds by the Markov inequality, the second inequality holds
by Assumption 4 and the fact that $\delta_{i}\leq 1$ for $i=1,\ldots,n$, the
third inequality holds since $\bm{X}_{1},\ldots,\bm{X}_{n}$ are identically
distributed, and the fourth inequality holds by Lemma 3. By (22), we have
validated (21), so we have

$\frac{1}{n}\sum_{i=1}^{n}\left\\{1-\frac{\delta_{i}}{\widehat{\pi}_{i}}\right\\}\left\\{\widehat{f}_{0}(\bm{x}_{i})-{f}^{*}(\bm{x}_{i})\right\\}=o_{\mathbb{P}}(n^{-1/2}).$
(23)

By (18), (19) and (23), we have proved Theorem 2 by the Slutsky’s theorem
(Athreya and Lahiri, 2006, Theorem 9.1.6). ∎

## References

  * (1)
  * Abadie and Imbens (2016) Abadie, A. and Imbens, G. W. (2016).  Matching on the estimated propensity score, Econometrica 84(2): 781–807.
  * Athreya and Lahiri (2006) Athreya, K. B. and Lahiri, S. N. (2006).  Measure Theory and Probability Theory, Springer, New York.
  * Banerjee et al. (2008) Banerjee, O., El Ghaoui, L. and d’Aspremont, A. (2008).  Model selection through sparse maximum likelihood estimation for multivariate Gaussian or binary data, Journal of Machine Learning Research 9(15): 485–516.
  * Barber and Candes (2019) Barber, R. and Candes, E. (2019).  A knockoff filter for high-dimensional selective inference, Annals of Statistics 47(5): 2504–2537.
  * Belloni et al. (2013) Belloni, A., Chernozhukov, V. and Hansen, C. (2013).  Inference on treatment effects after selection among high-dimensional controls, The Review of Economic Studies 81(2): 608–650.
  * Bertsimas and King (2017) Bertsimas, D. and King, A. (2017).  Logistic regression: From art to science, Statistical Science 32(3): 367–384.
  * Chen et al. (1999) Chen, M., Ibrahim, J. and Yiannoutsos, C. (1999).  Prior elicitation, variable selection and bayesian computation for logistic regression models, Journal of the Royal Statistical Society: Series B 61(1): 223–242.
  * Chernozhukov et al. (2018) Chernozhukov, V., Chetverikov, D., Demirer, M., Duflo, E., Hansen, C., Newey, W. and Robins, J. (2018).  Double/debiased machine learning for treatment and structural parameters, The Econometrics Journal 21(1): C1–C68.
  * Fan and Li (2001) Fan, J. and Li, R. (2001).  Variable selection via nonconcave penalized likelihood and its oracle properties, Journal of the American Statistical Association 96(456): 1348–1360.
  * Fan and Lv (2008) Fan, J. and Lv, J. (2008).  Sure independence screening for ultrahigh dimensional feature space (with discussion), Journal of the Royal Statistical Society, Series B 70(5): 849–911.
  * Fan et al. (2014) Fan, J., Xue, L. and Zou, H. (2014).  Strong oracle optimality of folded concave penalized estimation, Annals of statistics 42(3): 819.
  * Fan et al. (2011) Fan, J., Yang, F. and Song, R. (2011).  Nonparametric independence screening in sparse ultrahigh dimensional additive models, Journal of the American Statistical Association 106(494): 544–557.
  * Fischer and Steinwart (2020) Fischer, S. and Steinwart, I. (2020).  Sobolev norm learning rates for regularized least-squares algorithms, Journal of Machine Learning Research 21(205): 1–38.
  * Fukumizu and Leng (2014) Fukumizu, K. and Leng, C. (2014).  Gradient-based kernel dimension reduction for regression, Journal of the American Statistical Association 109(505): 359–370.
  * Genkin et al. (2007) Genkin, A., Lewis, D. D. and Madigan, D. (2007).  Large-scale Bayesian logistic regression for text categorization, Technometrics 49(3): 291–304.
  * Graham et al. (2012) Graham, B. S., de Xavier Pinto, C. C. and Egel, D. (2012).  Inverse probability tilting for moment condition models with missing data, The Review of Economic Studies 79(3): 1053–1079.
  * He et al. (2021) He, X., Wang, J. and Lv, S. (2021).  Efficient kernel-based variable selection with sparsistency, Statistica Sinica 31: 2123–2151.
  * Huang et al. (2010) Huang, J., Horowitz, J. and Wei, F. (2010).  Variable selection in nonparametric additive models, Annals of Statistics 38(4): 2282–2313.
  * Jaakkola et al. (1999) Jaakkola, T., Diekhans, M. and Haussler, D. (1999).  Using the Fisher kernel method to detect remote protein homologies, In Proceedings of Seventh International Conference on Intelligent Systems for Molecular Biology 99: 149–158.
  * Kim and Shao (2013) Kim, J. and Shao, J. (2013).  Statistical Methods for Handling Incomplete Data, CRC press, New York.
  * Kim et al. (2006) Kim, Y., Kim, J. and Kim, Y. (2006).  Blockwise sparse regression, Statistica Sinica 16(2): 375–390.
  * Lin et al. (2018) Lin, H., Zhou, F., Wang, Q., Zhou, L. and Qin, J. (2018).  Robust and efficient estimation for the treatment effect in causal inference and missing data problems, Journal of Econometrics 205: 363–380.
  * Lin and Zhang (2006) Lin, Y. and Zhang, H. (2006).  Component selection and smoothing in multivariate nonparametric regression, Annals of Statistics 34(5): 2272–2297.
  * Lv et al. (2018) Lv, S., Lin, H., Lian, H. and Huang, J. (2018).  Oracle inequalities for sparse additive quantile regression in reproducing kernel Hilbert space, Annals of Statistics 46(2): 781–813.
  * Meier et al. (2008) Meier, L., Van de Geer, S. and B$\ddot{{u}}$hlmann, P. (2008).  The group lasso for logistic regression, Journal of the Royal Statistical Society, Series B 70(1): 53–71.
  * Mendelson and Neeman (2010) Mendelson, S. and Neeman, J. (2010).  Regularization in kernel learning, Annals of Statistics 38(1): 526–565.
  * Mercer (1909) Mercer, J. (1909).  Functions of positive and negative type and their connection with the theory of integral equations, Philosophical Transactions of the Royal Society A 209: 415–446.
  * Ning et al. (2020) Ning, Y., Sida, P. and Imai, K. (2020).  Robust estimation of causal effects via a high-dimensional covariate balancing propensity score, Biometrika 107(3): 533–554.
  * Qin et al. (2017) Qin, J., Zhang, B. and Leung, D. H. (2017).  Efficient augmented inverse probability weighted estimation in missing data problems, Journal of Business & Economic Statistics 35(1): 86–97.
  * Qu et al. (2010) Qu, A., Lindsay, B. G. and Lu, L. (2010).  Highly efficient aggregate unbiased estimating functions approach for correlated data with missing at random, Journal of the American Statistical Association 105(489): 194–204.
  * Ravikumar et al. (2010) Ravikumar, P., Wainwright, M. and Lafferty, J. (2010).  High-dimensional Ising model selection using $\ell^{1}$-regularized logistic regression, Annals of Statistics 38(3): 1287 – 1319.
  * Robins et al. (1994) Robins, J. M., Rotnitzky, A. and Zhao, L. P. (1994).  Estimation of regression coefficients when some regressors are not always observed, Journal of the American Statistical Association 89(427): 846–866.
  * Rosasco et al. (2013) Rosasco, L., Villa, S., Mosci, S., Santoro, M. and Verri, A. (2013).  Nonparametric sparsity and regularization, Journal of Machine Learning Research 14(16): 1665–1714.
  * Rosenbaum and Rubin (1983) Rosenbaum, P. R. and Rubin, D. B. (1983).  The central role of the propensity score in observational studies for causal effects, Biometrika 70(1): 41–55.
  * Rubin (1976) Rubin, D. (1976).  Inference and missing data, Biometrika 63(3): 581–592.
  * Shen et al. (2012) Shen, X., Pan, W. and Zhu, Y. (2012).  Likelihood-based selection and sharp parameter estimation, Journal of the American Statistical Association 107(497): 223–232.
  * Shen et al. (2013) Shen, X., Pan, W., Zhu, Y. and Zhou, H. (2013).  On constrained and regularized high-dimensional regression, Annals of the Institute of Statistical Mathematics 65(5): 807–832.
  * Shevade and Keerthi (2003) Shevade, S. and Keerthi, S. (2003).  A simple and efficient algorithm for gene selection using sparse logistic regression, Bioinformatics 19(17): 2246–2253.
  * Smale and Zhou (2007) Smale, S. and Zhou, D. (2007).  Learning theory estimates via integral operators and their approximations, Constructive Approximation 26(2): 153–172.
  * Steinwart (2005) Steinwart, I. (2005).  Consistency of support vector machines and other regularized kernel classifiers, IEEE Transactions on Information Theory 51(1): 128–142.
  * Steinwart and Christmann (2008) Steinwart, I. and Christmann, A. (2008).  Support Vector Machine, Springer, New York.
  * Sun et al. (2013) Sun, W., Wang, J. and Fang, Y. (2013).  Consistent selection of tuning parameters via variable selection stability, Journal of Machine Learning Research 14(71): 3419–3440.
  * Tan (2010) Tan, Z. (2010).  Bounded, efficient and doubly robust estimation with inverse weighting, Biometrika 97(3): 661–682.
  * Tan (2020) Tan, Z. (2020).  Model-assisted inference for treatment effects using regularized calibrated estimation with high-dimensional data, Annals of Statistics 48(2): 811 – 837.
  * Tang et al. (2014) Tang, J., Alelyani, S. and Liu, H. (2014).  Feature selection for classification: A review, Data classification: Algorithms and applications pp. 37–64.
  * Tibshirani (1996) Tibshirani, R. (1996).  Regression shrinkage and selection via the lasso, Journal of the Royal Statistical Society, Series B 58(1): 267–288.
  * Van der Vaart (2000) Van der Vaart, A. (2000).  Asymptotic Statistics, Cambridge University Press, New York.
  * Wainwright et al. (2007) Wainwright, M., Ravikumar, P. and Lafferty, J. (2007).  High-dimensional graphical model selection using $\ell^{1}$-regularized logistic regression, Advances in neural information processing systems 19: 1465–1472.
  * Wang (2009) Wang, H. (2009).  Forward regression for ultra-high dimensional variable screening, Journal of the American Statistical Association 104(488): 1512–1524.
  * Wang and Kim (2021) Wang, H. and Kim, J. (2021).  Statistical inference after kernel ridge regression imputation under item nonresponse, Technical Report (Available at https://arxiv.org/pdf/2102.00058.pdf) pp. 1–26.
  * Wang and Leng (2016) Wang, X. and Leng, C. (2016).  High dimensional ordinary least squares projection for screening variables, Journal of the Royal Statistical Society, Series B 78(3): 589–611.
  * Wooldridge (2007) Wooldridge, J. (2007).  Inverse probability weighted estimation for general missing data problems, Journal of Econometrics 141(2): 1281–1301.
  * Yang et al. (2016) Yang, L., Lv, S. and Wang, J. (2016).  Model-free variable selection in reproducing kernel Hilbert space, Journal of Machine Learning Research 17(82): 1–24.
  * Yang et al. (2020) Yang, S., Kim, J. and Song, R. (2020).  Doubly robust inference when combining probability and non-probability samples with high-dimensional data, Journal of the Royal Statistical Society: Series B 82(2): 445–465.
  * Yuan and Lin (2006) Yuan, M. and Lin, Y. (2006).  Model selection and estimation in regression with grouped variables, Journal of the Royal Statistical Society: Series B 68(1): 49–67.
  * Zhang et al. (2016) Zhang, C., Liu, Y. and Wu, Y. (2016).  On quantile regression in reproducing kernel Hilbert spaces with data sparsity constraint, Journal of Machine Learning Research 17(40): 1–45.
  * Zhao et al. (2017) Zhao, P., Tang, N., Qu, A. and Jiang, D. (2017).  Semiparametric estimating equations inference with nonignorable missing data, Statistica Sinica 27(1): 89–113.
  * Zhou (2007) Zhou, D. (2007).  Derivative reproducing properties for kernel methods in learning theory, Journal of Computational and Applied Mathematics 220(1–2): 456–463.
  * Zou (2006) Zou, H. (2006).  The adaptive lasso and its oracle properties, Journal of the American Statistical Association 101(476): 1418–1429.