# Non-Asymptotic Guarantees for Robust Statistical Learning under
$(1+\varepsilon)$-th Moment Assumption

Lihu Xu<EMAIL_ADDRESS>
1\. Department of Mathematics, Faculty of Science and Technology, University
of Macau, Taipa Macau, China; 2\. Zhuhai UM Science & Technology Research
Institute, Zhuhai, China. Fang Yao<EMAIL_ADDRESS>
School of Mathematical Sciences and Center for Statistical Science, Peking
University, Beijing, China. Qiuran Yao<EMAIL_ADDRESS>
1\. Department of Mathematics, Faculty of Science and Technology, University
of Macau, Taipa Macau, China; 2\. Zhuhai UM Science & Technology Research
Institute, Zhuhai, China. Huiming Zhang<EMAIL_ADDRESS>
1\. Department of Mathematics, Faculty of Science and Technology, University
of Macau, Taipa Macau, China; 2\. Zhuhai UM Science & Technology Research
Institute, Zhuhai, China Corresponding authors: Huiming Zhang
(huimingzhang@um.edu.mo). The authors are in alphabetical order.

###### Abstract

There has been a surge of interest in developing robust estimators for models
with heavy-tailed data in statistics and machine learning. This paper proposes
a log-truncated M-estimator for a large family of statistical regressions and
establishes its excess risk bound under the condition that the data have
$(1+\varepsilon)$-th moment with $\varepsilon\in(0,1]$. With an additional
assumption on the associated risk function, we obtain an $\ell_{2}$-error
bound for the estimation. Our theorems are applied to establish robust
M-estimators for concrete regressions. Besides convex regressions such as
quantile regression and generalized linear models, many non-convex regressions
can also be fit into our theorems, we focus on robust deep neural network
regressions, which can be solved by the stochastic gradient descent
algorithms. Simulations and real data analysis demonstrate the superiority of
log-truncated estimations over standard estimations.

Keywords: data with $(1+\varepsilon)$-th moment, excess risk bounds, robust
quantile regression, robust generalized linear models (GLMs), robust non-
convex regressions, robust deep neural network (DNN) regressions.

##  1 Introduction

###  1.1 Backgrounds

Robust statistics is a traditional topic that has been well studied since the
pioneering work of Huber (1964) and Tukey (1960) in 1960s. Distributionally
robust learning is nowadays revitalized and invigorating in statistical
learning, see Nemirovskij and Yudin (1983) for median-of-means estimators,
Baraud et al. (2017) for minimax types estimators and Catoni (2012) for log-
truncated estimators. For further details, we refer the reader to the note by
Lerasle (2019) and the review paper by Lugosi and Mendelson (2019) for
comprehensive introductions.

Many existing pieces of research on excess risk bounds for robust estimators
heavily rely on one or more of the following three assumptions: (1) bounded
loss functions (Bartlett and Mendelson, 2006; Yi et al., 2020); (2) bounded
Lipschitz condition (Chinot et al., 2019; Shen et al., 2021a); (3) bounded
data (Liu and Tao, 2014; Brownlees et al., 2015; Zhang and Zhou, 2018) or
unbounded data with sub-Gaussian assumption (Xu et al., 2020; Lecué and
Mendelson, 2013; Loh, 2017; Ostrovskii and Bach, 2021). However, there are a
lot of statistical models which do not satisfy any of the above three
assumptions, see for instance Zhang and Jia (2022) and Chi (2010).

When a distribution does not have an exponential moment, it is often called
heavy-tailed; see Resnick (2007). However, there are many kinds of data, such
as network and wealth distribution data, which only have finite $p$-th moment
with $p\in(1,\infty)$. Zhang and Zhou (2018) studied log-truncated M-estimator
for least absolute deviation (LAD) regression under the assumption that the
data have $2$nd moment, while Chen et al. (2021a) extended their work to the
data with $(1+\varepsilon)$-th moment for $\varepsilon\in(0,1)$.

###  1.2 Contributions

This paper proposes a log-truncated M-estimator for a large family of
statistical regressions and establishes its excess risk bound under the
condition that the data have $(1+\varepsilon)$-th moment with
$\varepsilon\in(0,1]$, our _contributions_ are summarized as the following
three aspects.

  * •

The general high order function $\lambda(x)$ in (4) below. Zhang and Zhou
(2018) borrowed the idea from Catoni (2012) and choose $\lambda(x)=x^{2}/2$ to
study robust LAD regression when the data have $2$nd moment, while Chen et al.
(2021a) chose $\lambda(x)=\frac{|x|^{1+\varepsilon}}{1+\varepsilon}$ and
extended their results to the case in which the data have finite
$(1+\varepsilon)$-th moment. Table 1 below lists the choices of $\lambda(x)$
which has been reported in literatures. In Theorem 1, for finite
$(1+\varepsilon)$-th moment data with $\varepsilon\in(0,1]$, we find general
conditions for high order function $\lambda(x)$ (see (C.0)-(C.4) below) which
guarantee that the minimization of log-truncated loss function leads to a
consistent estimator in the sense of excess risk. With an additional
assumption on the risk function, we further show that the estimator is
consistent in $\ell_{2}$ norm.

  * •

General loss functions. Taking
$\lambda(x)=\frac{|x|^{1+\varepsilon}}{1+\varepsilon}$ with
$\varepsilon\in(0,1]$ and applying Theorem 1 with general loss functions, we
establish a new log-truncated robust estimator in Theorem 2, which not only
extends the results about LAD estimators in Zhang and Zhou (2018) and Chen et
al. (2021a), but also covers many other convex loss examples such as robust
quantile regression (QR) and robust generalized linear models (GLMs). For
GLMs, we obtain a general result for the related excess risk and apply it to
two typical models, logistic regression and negative binomial regression.

  * •

Non-convex problems: deep neural network regressions and beyond. Besides the
convex regressions mentioned above, we find that Theorem 2 can also be applied
to study non-convex regressions. This paper focuses on deep neural network
(DNN) regressions that can be solved by stochastic gradient descents (SGDs) in
practice and establish a general theorem. Empirical studies including Boston
housing, MNIST, and German health care demand datasets are performed well by
the proposed robust DNN regression models. We stress that Theorem 2 can be
applied to many other non-convex regressions, e.g. robust two-component mixed
linear regression and robust non-negative matrix factorization, in which one
has to design specific algorithms rather than use SGDs.

###  1.3 Related works

Catoni (2012) put forward a logarithm truncation for mean regression with
finite $2$nd moment data and obtained an estimator whose confidence interval
has a length comparable with that of the classical mean estimation with sub-
Gaussian data. Since then, Catoni’s idea has been extensively applied to study
regressions and estimations with heavy-tailed data; see Fan et al. (2017); Sun
et al. (2020); Wang et al. (2022); Sun (2021). Zhang and Zhou (2018) used
Catoni’s truncation technique to study a LAD regression for the data with
$2$nd moment and showed that the associated estimator is consistent in the
measure of excess risk. By modifying the logarithm truncation of Catoni, Chen
et al. (2021a) extended the results in Zhang and Zhou (2018) to the data with
$(1+\varepsilon)$-th moment with $\varepsilon\in(0,1)$. Due to the increasing
applications whose data do not have $2$nd moment, more and more
generalizations of Catoni’s truncation have been proposed, see Lam and Cheng
(2021), and Lee et al. (2020). Xu et al. (2020) studied the excess risk bounds
for learning with general non-convex truncated losses, in which
$\lambda(x)=O(x^{1+\varepsilon})$ with $\varepsilon=0$ or $\varepsilon\geq 1$.
Under bounded input assumption, Shen et al. (2021b) studied non-asymptotic
error bounds for deep neural networks for regression models with heavy-tailed
error output having a finite $(1+\varepsilon)$-th moment. For canonical GLMs,
Zhu and Zhou (2021) required $4$th moment condition on output to study the
consistency property of their proposed robust estimators.

###  1.4 Notations and organizations

The following notations will be frequently used in the rest of this paper.
Define the index set $[n]:=\\{1,2,\cdots,n\\}$. The r.v. is the shorthand for
random variable. Let $q\geq 1$, for $\theta\in\mathbb{R}^{p}$, define
${\|\theta\|_{q}}:={({\sum_{j=1}^{p}{\theta_{j}^{q}}})^{1/q}}$. Let
$(\Theta,d)$ be a metric space, and ${K}$ be a subset of $\Theta$. Define the
Euclidean norm ball $B_{r}=\left\\{x\in\mathbb{R}^{p}:\|x\|_{{2}}\leq
r\right\\}$. For an $\varepsilon>0$,
$\mathcal{N}(K,\varepsilon)\subset\mathbb{R}^{p}$ is an _$\varepsilon$ -net_
of $K$ if for all ${x}\in K,$ there is a ${y}\in\mathcal{N}(K,\varepsilon),$
such that $\|{x}-{y}\|_{2}\leq\varepsilon$. The _covering number_
${N}(K,\varepsilon)$ is the smallest number of closed balls centered at $K$
with radii $\varepsilon$ whose union covers $K$. For any scalar r.v. $Z$,
$\|\cdot\|_{\psi_{2}}$ is the sub-Gaussian norm defined by
$\|Z\|_{\psi_{2}}=\sup_{p\geq 1}\frac{({\rm{E}}|Z{|^{p}})^{1/p}}{\sqrt{p}}$.
For two matrices $A$ and $B$ with compatible dimension, denote $A\succ B$ if
$B-A$ is positive definite. Let
$\|A\|_{F}:=\sqrt{\sum_{i=1}^{m}\sum_{j=1}^{n}a_{ij}^{2}}$ be the Frobenius
norm of a square matrix $A=(a_{ij})\in\mathbb{R}^{m\times n}$.

The rest of the paper is organized as follows. Section 2 introduces the robust
estimator based on the log-truncated loss function and provides two general
theorems: Theorem 1 and Theorem 2, while Section 3 gives examples in which
several convex and non-convex regressions can be fit in the framework of
Theorem 2. Section 4 includes simulations and real data analysis, which
evaluate the effectiveness of the proposed log-truncated estimation for some
regressions discussed in Section 3, and Section 5 provides some perspectives
for the future study.

##  2 Log-truncated losses for data with $(1+\varepsilon)$-th moment

###  2.1 Problem setup

We assume that $\\{(X_{i},Y_{i})\\}_{i=1}^{n}$ are a sequence of
$\mathbb{R}^{p}\times\mathbb{R}$-valued independent identically distributed
(i.i.d.) r.v.s and that each $(X_{i},Y_{i})$ is a copy of r.v. $(X,Y)$. Denote
loss function by ${l(y,x,\theta)}$, where $y\in\mathbb{R}$ is the output
variable, $x\in\mathbb{R}^{p}$ is the input variable, and $\theta\in{\Theta}$
with $\Theta$ being the hypothesis space. Define

${{R}_{l}}(\theta):=\mathrm{E}[{l({Y},X,\theta)}],$

the true parameter $\theta^{*}\in\mathbb{R}^{p}$ is defined as the minimizer
of $R_{l}(\cdot)$, i.e.,

$\theta^{*}:=\mathop{\arg\min}\limits_{\theta\in{\Theta}}{{R}_{l}}(\theta),$
(1)

if $R_{l}(\theta^{*})<\infty$. The $\theta^{*}$ is called true parameter, it
may be not unique.

The estimator of $\theta^{*}$ based on the i.i.d. data
$\\{(X_{i},Y_{i})\\}_{i=1}^{n}$ is obtained by the following _empirical risk
minimization_ (ERM) problem:

$\bar{\theta}_{n}:=\mathop{\arg\min}\limits_{\theta\in{\Theta}}{{\hat{R}}_{l}^{n}}(\theta)\
\ \ \ \ \ \ {\rm with}\ \
{{\hat{R}}_{l}^{n}}(\theta)=\frac{1}{n}\sum\limits_{i=1}^{n}{l({Y_{i}},X_{i},\theta)}.$
(2)

Note that our loss function $l(\cdot,\cdot,\theta)$ may be non-convex, and
thus has many local minimizers in the above ERM problem. Various ERM examples
of non-convex loss appears in _Gaussian mixture models_ (Mei et al., 2018),
_mixture density estimation_ (Khamaru and Wainwright, 2019), _the mixture of
two linear regressions_ (Klusowski et al., 2019), _truncated Cauchy non-
negative matrix factorization_ (Guan et al., 2017), statistical learning under
deep neural networks (Fan et al., 2021).

We assume that our parameter space $\Theta$ is a bounded set in
$\mathbb{R}^{p}$, more precisely,

  * •

(C.0): $\Theta\subseteq\mathbb{R}^{p}$ is convex and there exists some
$r\in(0,\infty)$ such that $\|\theta\|_{2}\leq r,~{}\forall\theta\in\Theta$.

When the data are heavy tailed, the estimator $\bar{\theta}_{n}$ in (2) may
not be a robust estimator. Assuming that the loss function
$l(\cdot,\cdot,\theta)$ is _three times continuously differentiable self-
concordant_ and that the data are i.i.d. _sub-Gaussian_ , Ostrovskii and Bach
(2021) obtained an upper bound of the excess risk with high probability
$1-\delta$ , which has an order of
$\frac{O(1)}{{n}}\left({{{{\left\|{{\left\|{{X_{1}}}\right\|_{2}}}\right\|}_{{\psi_{2}}}}}}\log(\frac{e}{\delta})\right)$.
Obviously this bound will blow up when $X_{1}$ is not sub-Gaussian. Using
smoothed convex Huber losses and under the finite 8-$th$ moment assumption on
the data, Mathieu and Minsker (2021) studied regressions with quadratic losses
and obtained an excess risk bound with an order of
$O(\frac{p}{{n}}+\frac{1}{{n^{3/4}}})$. When the data do not have 8-$th$
moment, the method therein will not work.

We pursue other learning procedures by modifying the ERM to solve these
problems. In the sequel, we focus on log-truncated loss under the condition
that the data have $(1+\varepsilon)$-th moment.

###  2.2 Log-truncated loss and our estimator

Catoni (2012) put forward a non-deceasing truncation function $\psi$ such that

$-\log\left(1-x+\frac{x^{2}}{2}\right)\leq\psi(x)\leq\log\left(1+x+\frac{x^{2}}{2}\right),$
(3)

and obtained a robust mean estimator for i.i.d. data with finite second
moment. The truncation function $\psi(x)$ not only reduces the value of the
exponential-scaled outliers but also largely retains the data fluctuation in
an unbounded way, whereas the classical bounded truncated M-functions often
lose the information of the data with large values.

Using Catoni’s idea, (Zhang and Zhou, 2018) studied an LAD regression problem
under the same assumption. Slightly modifying the truncation function $\psi$,
Chen et al. (2021a) extended the results above to the data which have
$(1+\varepsilon)$-th moment with $\varepsilon\in(0,1)$. Minsker (2018), Lam
and Cheng (2021) and Lee et al. (2020) also considered the extensions of
$\psi$ for different motivations (see the table below), the key point for
these extensions is modifying the high order function $x^{2}/2$ in (3)
according to concrete problems at hand.

For $M$-estimation problems, we propose the following conditions for the high
order function $\lambda(x)$ such that the associated minimization of log-
truncated loss function leads to a consistent estimator in the sense of excess
risk. Let $\psi_{\lambda}:\mathbb{R}\rightarrow\mathbb{R}_{+}$ be a continuous
non-decreasing function such that

$\displaystyle-\log\left[1-x+\lambda(|x|)\right]\leq{\psi_{\lambda}}(x)\leq\log\left[1+x+\lambda(|x|)\right],\quad\forall
x\in\mathbb{R}$ (4)

where

  * •

(C.1): High order function
$\lambda(x):\mathbb{R}^{+}\rightarrow\mathbb{R}^{+}$ is a continuous non-
decreasing function such that
$\lim_{x\rightarrow\infty}\frac{\lambda(x)}{x}=\infty$. Moreover, there exist
some $c_{2}>0$ and a function $f:\mathbb{R}_{+}\rightarrow\mathbb{R}_{+}$ such
that

    * –

(C.1.1) $\lambda(tx)\leq f(t)\lambda(x)$ for all $t,x\in\mathbb{R}^{+}$, where
$f(t)$ is non-decreasing for $t\in(0,\infty)$ and $\lim_{t\rightarrow
0^{+}}{f(t)}/{t}=0$;

    * –

(C.1.2) $\lambda(x+y)\leq c_{2}[\lambda(x)+\lambda(y)]$ for all
$x,y\in\mathbb{R}^{+}$.

We call the inequality in (C.1.1) as _homogeneous inequality_ and the one in
(C.1.2) as _weak triangle inequality_.

Inspired by (Zhang and Zhou, 2018; Chen et al., 2021a), we define an estimator
$\hat{\theta}_{n}$ by the following ERM problem:

$\hat{\theta}_{n}:=\mathop{\arg\min}\limits_{\theta\in{\Theta}}{\hat{R}}^{n}_{\psi_{\lambda},l,\alpha}(\theta)\
\ \ {\rm with}\ \ \ \ \ \
{\hat{R}}^{n}_{\psi_{\lambda},l,\alpha}(\theta):=\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}{\psi_{\lambda}}[\alpha
l({Y_{i}},{{X}}_{i},\theta)],$ (5)

where $\alpha>0$ is a parameter to be tuned. The problem (5) is non-convex
since the $\psi_{\lambda}$ is not convex. For further use, we define the
following risk function:

${R_{\lambda\circ l}}({\theta}):={{\rm{E}}}\\{\lambda[l({Y},X,\theta)]\\},\ \
\ \ \theta\in\Theta.$

Table 1: Some examples of high order function $\lambda(|x|)$ and $\varepsilon\in(0,1)$. Higher order functions $\lambda(|x|)$ | References | Moment conditions
---|---|---
$\frac{1}{2}|x|^{2}$ | Catoni (2012) | $\mathrm{E}\left|X\right|^{2}<\infty$
$\left(\frac{\varepsilon}{1+\varepsilon}\vee\sqrt{\frac{1-\varepsilon}{1+\varepsilon}}\right)|x|^{1+\varepsilon}$ | Minsker (2018) | $\mathrm{E}\left|X\right|^{1+\varepsilon}<\infty$
$\frac{1}{1+\varepsilon}|x|^{1+\varepsilon}$ | Chen et al. (2021a) | $\mathrm{E}\left|X\right|^{1+\varepsilon}<\infty$
${\left[{2{{(\frac{{1-\varepsilon}}{\varepsilon})}^{1-2/(1+\varepsilon)}}+{{(\frac{{1-\varepsilon}}{\varepsilon})}^{2-2/(1+\varepsilon)}}}\right]^{-(1+\varepsilon)/2}}|x|^{1+\varepsilon}$ | Lee et al. (2020) | $\mathrm{E}\left|X\right|^{1+\varepsilon}<\infty$
$(1+\varepsilon)^{-(1+\varepsilon)/2}({1-\varepsilon})^{1-(1+\varepsilon)/2}\varepsilon^{\varepsilon}{|x|^{1+\varepsilon}}$ | Lam and Cheng (2021) | $\mathrm{E}\left|X\right|^{1+\varepsilon}<\infty$
$\sum_{k=2}^{m}\frac{x^{k}}{k!},~{}(m\geq 2)$ | Xu et al. (2020) | $\sum_{k=2}^{m}\frac{\mathrm{E}\left|X\right|^{k}}{k!}<\infty$

Let us give an intuition of choosing $\psi_{\lambda}$ and $\alpha$. Without
loss of generality, we choose $\psi_{\lambda}(x)$ in (4) as follows:

$\displaystyle{\psi_{\lambda}}(x)$
$\displaystyle=\log\left[{1+x+\lambda(|x|)}\right]{\rm{1(}}x\geq
0{\rm{)}}-\log\left[{1-x+\lambda(|x|)}\right]{\rm{1(}}x\leq 0{\rm{)}}$
$\displaystyle={\rm{sign}}(x)\log(1+|x|+\lambda(|x|)).$ (6)

For example, let
$\lambda(x)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon},~{}{\varepsilon}\in(0,1]$
and the score function is
$\dot{\psi}_{\lambda}(x)={\frac{{1+|x|^{\varepsilon}}}{{1+|x|+|x|^{1+\varepsilon}/(1+\varepsilon)}}}\in(0,1]$.
Thus

$\displaystyle\nabla{\hat{R}_{\psi_{\lambda},l,\alpha}}(\theta)=\frac{1}{n}\sum\limits_{i=1}^{n}{X_{i}\dot{l}({Y_{i}},X_{i}^{\top}\theta)\dot{\psi}}[\alpha
l({Y_{i}},X_{i}^{\top}\theta)]\approx\frac{1}{n}\sum\limits_{i=1}^{n}X_{i}\dot{l}({Y_{i}},X_{i}^{\top}\theta),$
(7)

where the approximation is by putting

$\displaystyle\alpha\to 0~{}\text{and thus}~{}{\dot{\psi}[\alpha
l({Y_{i}},X_{i}^{\top}{\theta^{*}})]}\to 1~{}\text{as}~{}n\to\infty.$ (8)

The $\alpha$ will be chosen according to the sample size $n$ and tend to $0$
as $n\rightarrow\infty$. In order to obtain the optimal ${\alpha}$ as deriving
excess risk bounds, we assume the following conditions for further analysis:

  * •

(C.2) Local Lipschitz condition: $\exists$ _locally Lipschitz constant_
$H(y,x)$ s.t. $l(y,x,\cdot)$ satisfies

$|l(y,x,{\eta_{2}})-l(y,x,{\eta_{1}})|\leq
H(y,x)\|{\eta_{2}}-{\eta_{1}}\|_{2};{\eta_{1}},{\eta_{2}}\in\Theta.$

  * •

(C.3) Moment condition: ${{\rm{E}}}\\{\lambda(H(Y,X))\\}<\infty$.

  * •

(C.4) The existence of risk functions:
${\sup}_{\theta\in\Theta}R_{\lambda\circ l}(\theta)<\infty$.

The (C.2) can be obtained by a first-order Taylor expansion of $l(y,x,\cdot)$
as the following

$l(y,x,{\eta_{2}})=l(y,x,{\eta_{1}})+({\eta_{2}}-{\eta_{1}})^{\top}\dot{l}[y,x,(t{\eta_{2}}+(1-t){\eta_{1}})];{\eta_{1}},{\eta_{2}}\in\Theta,~{}\exists~{}t\in(0,1),$

where $\dot{l}(y,x,\cdot)$ is a (sub-)gradient, we can choose a $H(y,x)$
satisfying

$H(y,x)\geq{\sup}_{{\eta_{1}},{\eta_{2}}\in\Theta}\|{\dot{l}[y,x,(t{\eta_{2}}+(1-t){\eta_{1}})]}\|_{2}.$
(9)

The (C.2) is also indispensable for controlling remainder terms in the proof
of asymptotic normality of M-estimators, by showing the uniformly convergence
of empirical processes concerning the M-functions; see Section 19.4 and Lemma
19.31 in Van der Vaart (1998) for details. The (C.4) ensures that the excess
risk is well-defined.

###  2.3 Main results for learning with log-truncated losses

We are at the position to give _the first main result_ , in which the $\alpha$
is a turning parameter to be determined in concrete examples.

###### Theorem 1

Let $\hat{\theta}_{n}$ be defined by (5) with high order function
$\lambda(x)$, and let ${\theta^{*}}$ be a true value given by (1). For
constants $\delta\in(0,1/2)$ and $\kappa>0$, under (C.0)-(C.4), we have with
probability at least $1-2\delta$

$\displaystyle~{}~{}~{}~{}{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq
2\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{{(c_{2}+1)f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}{R_{\lambda\circ
l}(\theta)}+\frac{1}{{n\alpha}}\log\frac{N(\Theta,\kappa)}{\delta^{2}},$

where $c_{2}$ is a constant in (C.1.1) and $f(t)$ is the function in (C.1.2).
In particular, let $\kappa=1/n$,
${\underline{\alpha}_{n}}:={f^{-1}}\left({\frac{{\log({\delta^{-2}})+p\log\left({rn}\right)}}{{n({c_{2}}+1){\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)}}}\right)$ and
${{\bar{\alpha}}_{n}}:={f^{-1}}\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}}{{n({c_{2}}+1){\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)}}}\right)$, then

$R_{l}(\hat{\theta}_{n})-R_{l}(\theta^{*})\leq\frac{{2{{\rm{E}}}[H(Y,X)]}}{{n}}+\frac{{{c_{2}}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}}}{{{\underline{\alpha}_{n}}}}f(\frac{{\bar{\alpha}}_{n}}{n})+\frac{2[{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}]}{{n{\underline{\alpha}_{n}}}},$
(10)

with probability at least $1-2\delta$.

The terms
$\frac{{(c_{2}+1)f(\alpha)}}{\alpha}{\sup}_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)$ and
$\frac{1}{{n\alpha}}\log\frac{{N({\Theta},\kappa)}}{{{\delta^{2}}}}$ in
Theorem 1 can be viewed as variance and bias respectively. There is a tradeoff
between the confidence probability $1-2\delta$ and the bias
$\frac{1}{{n\alpha}}\log\frac{{N({\Theta},\kappa)}}{{{\delta^{2}}}}$, i.e.,
the larger the $1-2\delta$ is, the bigger the bias will be. We apply the
variance-bias tradeoff principle, i.e.
$\frac{{(c_{2}+1)f(\alpha)}}{\alpha}{\sup}_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)=\frac{1}{{n\alpha}}\log\frac{{N({\Theta},\kappa)}}{{{\delta^{2}}}}$,
to choose the tuning parameter $\alpha$. Note that Zhang and Zhou (2018) and
Theorem 4.1 in Chen et al. (2021a) chose their tuning parameters without
considering this principle. For i.i.d. data, Zhang and Zhou (2018) studied LAD
regression with $\lambda(x)=|x|^{2}/2$, while Chen et al. (2021a) considered
LAD regression with $\lambda(x)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon}$
and $\varepsilon\in(0,1)$. The following theorem extends these two results to
a large class of loss functions, which cover quantile regression, generalized
linear models, and non-convex regressions such as DNN regressions.

###### Theorem 2

By setting
$\lambda(x)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon},~{}{\varepsilon}\in(0,1]$
and
$\alpha\leq\frac{1}{n^{1/({1+\varepsilon})}}(\frac{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}{[(2^{\varepsilon}+1)\sup_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]})^{1/({1+\varepsilon})}$ in Theorem 1, with probability at least
$1-2\delta,~{}(\delta\in(0,1/2))$ one has

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})\leq\frac{2}{n}{{\rm{E}}}[H({Y},{{X}})]+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}$
$\displaystyle\times\left[\frac{2({2^{\varepsilon}}+1)}{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{-{({1+\varepsilon})^{-1}}}}+\frac{{{2^{\varepsilon}}{{\rm{E}}}|{H}(Y,X)|^{1+\varepsilon}/(1+\varepsilon)}}{{{n^{\varepsilon+1}}{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right]=O\left({{{\left({\frac{{p\log
n}}{n}}\right)}^{{\varepsilon}/{({1+\varepsilon})}}}}\right),$

where $p=p_{n}$ is allowed to grow with $n$. Additionally, if $(C.3)$ is
replaced by

${{\rm{E}}}[{H}(Y,X)]=o(n),{{\rm{E}}}|{H}(Y,X)|^{1+\varepsilon}=O({n^{1+\varepsilon}}).$

Then, the consistency of excess risk is valid under $p=o(n/\log n)$, i.e.

${R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})=o_{p}(1).$

The finite-sample upper bound for
${R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$ in Theorem 2 tells us that
the convergence rate of excess risk is of an order $O({{{({\frac{{p\log
n}}{n}})}^{{\varepsilon}/{({1+\varepsilon})}}}})$. The smaller $\varepsilon$
is, the slower convergence rate will be. In Section 3, we will illustrate
several applications of Theorem 2 under the finite $(1+\varepsilon)$-th moment
with $\varepsilon\in(0,1]$, including quantile regressions, GLMs, and robust
non-convex regression such as robust DNN regressions.

By definition of the true risk minimization (1), we have
${R_{l}}(\theta)-{R_{l}}({\theta^{*}})\geq 0$ for all $\theta\in\Theta$.
Further, we assume that the excess function satisfies the following lower
bound:

  * •

(C.5) Hessian condition for the risk function: $\exists~{}c_{l}>0$,

${R_{l}}(\theta)-{R_{l}}({\theta^{*}})\geq
c_{l}{\left\|{\theta-\theta^{*}}\right\|_{2}^{2}}$ for all $\theta\in\Theta$.

In general M-estimators, the (C.5) is used in (M1) of Shao and Zhang (2021)
for deriving Berry-Esseen bound, and it is also considered in Assumption 3.3
of Schreuder et al. (2020) for obtaining non-asymptotic LIL. In fact, (C.5) is
motivated from the following intuition. Given $\theta\in\Theta$, applying
Taylor expansion to ${R_{l}}(\theta)$ at ${\theta^{*}}$ gives

$\displaystyle{R_{l}}(\theta)-{R_{l}}(\theta^{*})={(\theta-\theta^{*})^{\top}}{\nabla^{2}{R}_{l}}(t\theta^{*}+(1-t)\theta)(\theta-\theta^{*}),$

where $t\in[0,1]$ and the second equality is by $\nabla{R}_{l}(\theta^{*})=0$.
Under (C.0), the (C.5) is naturally valid when the Hessian matrix
${\nabla^{2}{R}_{l}}(t\theta^{*}+(1-t)\theta)\succ c_{l}\mathrm{I}_{p}$ for
$t\in[0,1]$ under the restriction
$\left\|{t\theta^{*}+(1-t)\theta}\right\|_{2}\leq
t\left\|{\theta^{*}}\right\|_{2}+(1-t)\left\|{\theta}\right\|_{2}\leq r$.

Combing with (C.5), Theorem 1 immediately implies the $\ell_{2}$-risk bounds,
which guarantees the consistency of the proposed robust estimators.

######  Corollary 3 ($\ell_{2}$-estimation error)

Under the same condition in Theorem 1, if (C.5) holds additionally, we have
with probability at least $1-2\delta,~{}(\delta\in(0,1/2))$,

${\|{\hat{\theta}_{n}-\theta^{*}}\|_{2}^{2}}\leq\frac{{2{{\rm{E}}}[H(Y,X)]}}{{n{c_{l}}}}+{{\rm{E}}}\\{\lambda[H(Y,X)]\\}\frac{{c_{2}}{f({{\bar{\alpha}}_{n}}/n)}}{{{\underline{\alpha}_{n}}{c_{l}}}}+\frac{2[{\log({\delta^{-2}})+p\log\left({1+2r/\kappa}\right)}]}{{n{\underline{\alpha}_{n}}{c_{l}}}}.$

In particular, if we consider the setting in Theorem 2, then
${\|{\hat{\theta}_{n}-\theta^{*}}\|_{2}^{2}}\leq O(C_{\delta}{{{({\frac{{p\log
n}}{n}})}^{{\varepsilon}/{({1+\varepsilon})}}}})$ with probability at least
$1-2\delta$, for a certain constant $C_{\delta}>0$.

##  3 Examples

This section provides two types of examples: robust convex regression and
robust non-convex regression. Convex regression includes quantile regression
and GLMs, while non-convex regressions contain deep neural network regression
for LAD, logistic, and negative binomial regression models. We assume that the
data in this section has the finite $(1+\varepsilon)$-th moment with
$\varepsilon\in(0,1]$.

###  3.1 Robust regressions with convex losses

####  3.1.1 Quantile regression for non-identically distributed noises

Quantile regression (QR) is a regression that concerns on estimating the
conditional quantiles of the response, which is put forward by Koenker and
Bassett Jr (1978) as an extension of linear regression. Consider

$Y_{i}={{X}}_{i}^{\top}\theta^{*}+\epsilon_{i},(i=1,\ldots,n),$ (11)

where ${{X}}_{i}=(X_{i1},\ldots,X_{ip})^{\top}$ is the $i$-th stochastic
design point in $\mathbb{R}^{p}$, and $\epsilon_{i}$’s is random errors with
different probability density $f_{i}(\cdot)$’s. Here, we may assume that
$\epsilon_{i}$’s are independent and satisfy
$P\left(\epsilon_{i}<0|X_{i}\right)=\tau$ for $0<\tau<1$. The unknown
regression coefficient
$\theta^{*}=(\theta^{*}_{1},\ldots,\theta^{*}_{p})^{\top}$ may depend on
$\tau$ and sample size $n$, but we omit such dependence in notation for
simplicity. The conditional distribution of $Y$ given $x$ is $F(y|x)=P(Y\leq
y|x)$ and the $\tau\mathrm{th}$ conditional quantile of $Y$ given $x$ is
$Q_{y|x}(\tau)=\inf\\{t:F(t|x)\geq\tau\\}$. Thus the problem of interest is to
estimate the unknown slope coefficient $\theta^{*}$ by regressing the
conditional quantile function

$Q_{Y_{i}|X_{i}}(\tau)={{X}}_{i}^{\top}\theta^{*},(i=1,\ldots,n)$.

Recall that the loss function of quantile regression is

$l(y,x,\theta)={\rho_{\tau}}(y-x^{\top}\theta)$ with
$\rho_{\tau}(u)=u[\tau-I(u<0)]$.

Under the i.i.d. data $\\{(X_{i},Y_{i})\\}_{i=1}^{n}$, we study the log-
truncated estimator $\hat{\theta}_{n}$ for the quantile loss:

$\hat{\theta}_{n}:=\mathop{\arg\min}\limits_{\theta\in{\Theta}}{{\hat{R}}_{\psi_{\lambda},{\rho_{\tau}},\alpha}^{n}}(\theta),$
(12)

where
${{\hat{R}}_{\psi_{\lambda},{\rho_{\tau}},\alpha}^{n}}(\theta):=\frac{1}{{n\alpha}}\sum_{i=1}^{n}\psi_{\lambda}[\alpha{\rho_{\tau}}({Y_{i}}-{{X}}_{i}^{\top}\theta)]$
and
$\lambda(x)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon},~{}{\varepsilon}\in(0,1]$.
The _tuning parameter_ $\alpha$ will be specified. The true parameter
$\theta^{*}$ is defined as the minimizer

$\theta^{*}:=\mathop{\arg\min}\limits_{\theta\in{\Theta}}{{R}_{\rho_{\tau}}}(\theta),$
(13)

where
${{R}_{\rho_{\tau}}}(\theta):=\mathrm{E}[{{\rho_{\tau}}({Y}-X^{\top}\theta)}]$
for ${\theta\in{\Theta}}$. Besides (C.0), we further assume

  * •

(Q.1): ${{\rm{E}}}\left\|{{X}}\right\|_{{2}}^{1+\varepsilon}<\infty$.

  * •

(Q.2): Let
${{R}_{\lambda\circ\rho_{\tau}}}(\theta):=\mathrm{E}[\lambda({\rho_{\tau}}({Y}-{{X}}^{\top}\theta))]$.
We have
$\sup_{\theta\in\Theta}{{R}_{\lambda\circ\rho_{\tau}}}(\theta)<\infty$.

The (Q.1) indicates that finite $(1+\varepsilon)$-th moment for Euclidian norm
of input is sufficient. The (Q.2) implies (C.4). To guarantee high probability
$\ell_{2}$-estimation error, we need the following two extra conditions:

  * •

(Q.3): The ${\varepsilon}_{1}$ is the continuous r.v. and its density function
are uniformly bounded away from 0:
$\mathop{\inf}\limits_{|t|<\infty}{f_{\varepsilon_{1}}}(t)>{c_{f}},$ where
${c_{f}}$ is a positive constant.

  * •

(Q.4): Assume that the matrix $\mathrm{E}[XX^{\top}]$ is non-singular. That is
to say there exists a positive constant ${c_{\Sigma}}$ such that
${\rm{E}}[XX^{\top}]\succ{c_{\Sigma}}\mathrm{I}_{p}.$

By the Taylor-like expansion (Knight’s identity, Knight (1998)) of the check
function

$\rho_{\tau}(u-v)-\rho_{\tau}(u)=-v[\tau-{\rm{1}}(u<0)]+\int_{0}^{v}[{\rm{1}}(u\leq
s)-{\rm{1}}(u\leq 0)]ds,~{}u,v\in\mathbb{R}$

from (46) and (A.7) in Appendix, one can check
$H(y,x)=l_{\tau}\left\|x\right\|_{2}$ under (C.0), where
$l_{\tau}:=\max\\{1+\tau,2-\tau\\}$. Thus $H(y,x)$ is free of $y$ and Theorem
2 implies

###### Corollary 4

Let $\tau\in(0,1)$, $\delta\in(0,1/2)$ and $\varepsilon>0$. Define
$\hat{\theta}_{n}$ by (12), and ${\theta^{*}}$ is given by (13). Under (C.0),
(Q.1) and (Q.2), if we put
$\alpha\leq\frac{1}{n^{1/(\varepsilon+1)}}\left(\frac{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}{({2^{\varepsilon}}+1)\sup_{\theta\in\Theta}R_{\lambda\circ\rho_{\tau}}^{n}(\theta)}\right)^{1/(\varepsilon+1)}$
and
$C_{\delta,n,\varepsilon}:={\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}$.
Then, with probability at least $1-2\delta$ one has

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})\leq\frac{{2{l_{\tau}}}}{n}{\rm{E}}{\left\|{{X}}\right\|_{{2}}}+\frac{C_{\delta,n,\varepsilon}}{{{n^{\varepsilon/(\varepsilon+1)}}}}\left[\frac{2({2^{\varepsilon}}+1)}{{[\sup\limits_{\theta\in\Theta}R_{\lambda\circ\rho_{\tau}}^{n}(\theta)]}^{-{(\varepsilon+1)^{-1}}}}+\frac{{{2^{\varepsilon}(l_{\tau}/n)^{1+\varepsilon}}{{\rm{E}}}{\left\|{{X}}\right\|_{{2}}^{1+\varepsilon}}}}{{{{[{\sup\limits_{\theta\in\Theta}}R_{\lambda\circ\rho_{\tau}}^{n}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right].$

If (Q.3) and (Q.4) hold additionally, we have with probability at least
$1-2\delta$

${\|{\hat{\theta}_{n}-\theta^{*}}\|_{2}^{2}}\leq\frac{{2{l_{\tau}}}}{{{c_{f}}{c_{\Sigma}}}n}{\rm{E}}{\left\|{{X}}\right\|_{{2}}}+\frac{C_{\delta,n,\varepsilon}}{{{{c_{f}}{c_{\Sigma}}}{n^{{\varepsilon/(\varepsilon+1)}}}}}\left[\frac{2({2^{\varepsilon}}+1)}{{[\sup\limits_{\theta\in\Theta}R_{\lambda\circ\rho_{\tau}}^{n}(\theta)]}^{-{(\varepsilon+1)^{-1}}}}+\frac{{{2^{\varepsilon}(l_{\tau}/n)^{1+\varepsilon}}{{\rm{E}}}{\left\|{{X}}\right\|_{{2}}^{1+\varepsilon}}}}{{{{[{\sup\limits_{\theta\in\Theta}}R_{\lambda\circ\rho_{\tau}}^{n}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right].$

It should be noted that Section 4.1.2 in Koenker (2005) mentioned that $2$nd
moment condition of the input is required to show the consistency for the ERM
estimator
$\bar{\theta}_{n}:={\arg\min}_{\theta\in{\Theta}}\frac{1}{{n}}\sum_{i=1}^{n}{\rho_{\tau}}({Y_{i}}-{{X}}_{i}^{\top}\theta).$
When the sample $\\{(X_{i},Y_{i})\\}_{i=1}^{n}$ satisfies that
$\\{X_{i}\\}_{i=1}^{n}$ are i.i.d. while $\\{Y_{i}\\}_{i=1}^{n}$ are non-
identically distributed obeying (11) with heterogeneous regression noises [see
He and Shao (1996) and Section 3.2.3 in Koenker (2005)], similar version of
Corollary 4 with second moment assumption is presented in Appendix.

####  3.1.2 Robust generalized linear models (GLMs)

Moving beyond the quantile loss for linear models, we focus on general loss
function of GLMs (McCullagh and Nelder, 1989) as below. In this part, we
assume that $\\{Y_{i}\\}_{i=1}^{n}$ satisfy some moment conditions rather than
the assumptions for distribution (as in classical GLMs).

Let $u(\cdot)$ be a known link function. Consider the quasi-GLMs loss
function:

${{\hat{R}}_{l}^{n}}(\theta):=\frac{1}{n}\sum_{i=1}^{n}l({Y_{i}},X_{i}^{\top}\theta),~{}\theta\in\mathbb{R}^{p}$
(14)

where the loss function is given by
$l(y,x^{\top}\theta):=k(x^{\top}\theta)-yu(x^{\top}\theta)$ with
$k(t):=b(u(t))$ for a given function $b(\cdot)$. If $u(t)=t$, we say the
quasi-GLMs has the _canonical link_.

Let $\alpha$ be a tuning parameter to be specified, the log-truncated robust
estimator $\hat{\theta}$ for quasi-GLMs is

$\hat{\theta}_{n}:=\mathop{\arg\min}\limits_{\theta\in{\Theta}}{{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\theta),$
(15)

where
${{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\theta):=\frac{1}{{n\alpha}}\sum_{i=1}^{n}\psi_{\lambda}[\alpha(k(X_{i}^{\top}\theta)-Y_{i}u(X_{i}^{\top}\theta))]$
with
$\psi_{\lambda}(x):={\rm{sign}}(x)\log(1+|x|+\frac{1}{1+\varepsilon}|x|^{1+\varepsilon})$.
We assume the GLMs-related conditions:

  * •

(G.1): Assume that $u(\cdot)$ is continuous differentiable and
$\dot{u}(\cdot)\geq 0$, there exist a positive constant $A$ and a positive
function ${g_{A}}(\cdot)$:

$0\leq\dot{u}(x^{\top}\theta)\leq{g_{A}}(x),~{}\text{for}~{}\|\theta\|_{2}\leq
A.$

  * •

(G.2): Given a function $b(\cdot)$ such that $\ddot{b}(\cdot)>0$, suppose that
$k(\cdot)$ is continuous differentiable and $k(\cdot)\geq 0$, there exist a
positive function ${h_{A}}(\cdot)$:

$0<\dot{k}(x^{\top}\theta)\leq{h_{A}}(x),~{}\text{for}~{}\|\theta\|_{2}\leq
A.$

  * •

(G.3): Mixed moment condition:
${\rm{E}}|[|Y|{g_{A}}(X)+{h_{A}}(X)]\left\|X\right\|_{2}|^{1+\varepsilon}<\infty$.

  * •

(G.4): Assume
$\sigma_{R}:={\sup}_{\theta\in{\Theta}}{{\rm{E}}}[k(X^{\top}\theta)-Yu(X^{\top}\theta)]^{1+\varepsilon}/{(1+\varepsilon)}<\infty$.

The conditions (G.1), (G.2) and (G.3) imply (C.3), while (G.4) implies (C.4);
see Remark 12 in Appendix for more discussion for (G.1) and (G.2). Theorem 2
is applicable to obtain the following result.

###### Corollary 5

Let ${\theta^{*}}$ be defined by (1) with loss
$l(y,x,\theta):=k(x^{\top}\theta)-yu(x^{\top}\theta)$ defined in (14), and
$\hat{\theta}_{n}$ is given by (15). Under (C.0), (G.1)-(G.4), if
$\alpha\leq\frac{1}{n^{1/({1+\varepsilon})}}(\frac{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}{(2^{\varepsilon}+1)\sigma_{R}})^{1/({1+\varepsilon})}$,
then with probability at least $1-2\delta$

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{2}{n}{{\rm{E}}}\\{[|Y|{g_{A}}(X)+{h_{A}}(X)]\left\|X\right\|_{2}\\}+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}$
$\displaystyle\times\left[{2({2^{\varepsilon}}+1){\sigma_{R}^{{({1+\varepsilon})^{-1}}}}}+\frac{{{2^{\varepsilon}}{{\rm{E}}}|[|Y|{g_{A}}(X)+{h_{A}}(X)]\left\|X\right\|_{2}|^{1+\varepsilon}}}{{{(1+\varepsilon)n^{\varepsilon+1}}{\sigma_{R}^{\varepsilon/(\varepsilon+1)}}}}\right].$

Corollary 5 can be applied to the following two examples, robust logistic
regression and robust negative binomial regression.

######  Example 1 (Robust logistic regression)

The output in logistic regression can take only two values: “0, 1”. Formally,
let $Y_{i}$’s $\in\\{0,1\\}$ be the random outputs and ${\theta^{\rm{*}}}$ be
a $p\times 1$ vector of unknown regression coefficients belonging to a compact
subset of $\mathbb{R}^{p}$. Given $n$ random input $X_{i}$’s
$\in\mathbb{R}^{n\times p}$, the logistic regression assumes
$P({Y_{i}}=1|{X_{i}};{\theta^{\rm{*}}}):=\frac{{{e^{X_{i}^{\top}{\theta^{\rm{*}}}}}}}{{1+{e^{X_{i}^{\top}{\theta^{\rm{*}}}}}}}.$
The empirical loss function of logistic regression is

$\hat{R}_{l}(\theta)=\frac{-1}{n}\sum_{i=1}^{n}{[{Y_{i}}X_{i}^{\top}{\theta}-\log(1+{e^{X_{i}^{\top}{\theta}}})}]$.

To obtain the finite excess risk, the robust logistic regression requires the
moment condition

${\rm{E}}\|{{X}}\|_{2}^{1+\varepsilon}<\infty.$

We obtain

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{{4}{{\rm{E}}}{\left\|{X}\right\|_{2}}}{n}+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}\cdot$
$\displaystyle\left[\frac{2({2^{\varepsilon}}+1)}{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{-{(\varepsilon+1)^{-1}}}}+\frac{{2^{1+2\varepsilon}}{{\rm{E}}}{\left\|{X}\right\|_{2}^{1+\varepsilon}/{({1+\varepsilon})}}}{{n^{1+\varepsilon}{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right].$

For modeling count data regressions, the classical Poisson regression as the
canonical link GLMs has equal dispersion assumption (i.e.
$\mathrm{E}(Y|X)=\mathrm{Var}(Y|X)$), which has little practical motivation.
Nevertheless, it motivates us to study the more flexible count data
regressions, as shown below.

######  Example 2 (Robust negative binomial regression)

As a generalization of Poisson regression, negative binomial regression (NBR)
relaxes the equadispersion assumption to the quadratic relationship between
the mean and variance of the response. NBR assumes that the overdispersed
responses $\\{Y_{i}\\}_{i=1}^{n}$ are modelled by two-parameter negative
binomial distribution with the connection of covariates:
$P({Y_{i}=y}|\theta,\mu_{i})=\frac{\Gamma(\eta+y)}{\Gamma(\eta)y!}(\frac{\mu_{i}}{\eta+\mu_{i}})^{y}(\frac{\eta}{\eta+\mu_{i}})^{\eta}~{}\text{with}~{}\log\mu_{i}=X_{i}^{\top}\theta,$
where $\eta>0$ is the known dispersion parameter, which can be estimated
previously. We have
${\rm{E}}({Y_{i}}|{X}_{i})={\mu_{i}}\leq{\rm{Var}}({Y_{i}}|X_{i})={\mu_{i}}+{{\mu_{i}^{2}}}/{\eta}$.
The NBR loss function is

${{\hat{R}}_{l}^{n}}(\theta)=\frac{-1}{n}\sum\limits_{i=1}^{n}\\{Y_{i}[{X}_{i}^{\top}\theta-\log(\eta+e^{{X}_{i}^{\top}\theta})]-\eta\log(\eta+e^{{X}_{i}^{\top}\theta})\\},$

see Zhang and Jia (2022) for more details. From Corollary 5,

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{2{{\rm{E}}}[(Y+\eta)\left\|X\right\|_{2}]}{n}+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}\cdot$
$\displaystyle\left[\frac{2({2^{\varepsilon}}+1)}{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{-{(\varepsilon+1)^{-1}}}}+\frac{{{{{4^{\varepsilon}}}{{\rm{E}}}[{\left\|YX\right\|_{2}^{1+\varepsilon}}+{(\eta\left\|X\right\|_{2})^{1+\varepsilon}/{({1+\varepsilon})}}]}}}{n^{1+\varepsilon}{{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right].$

Note that there are no assumptions for the distribution of output, and it only
requires the moment conditions:

${\rm{E}}\|{{X}{Y}}\|_{2}^{1+\varepsilon}<\infty$ and
${\rm{E}}\|{{X}}\|_{2}^{1+\varepsilon}<\infty$.

###  3.2 Robust regressions with non-convex losses: deep neural network

In many statistical learning problems, loss functions are non-convex due to
the high complexity of the models, and the ERM will have multiple local
minima. Classical examples of non-convex problems include matrix factorization
(Guan et al., 2017), principle component analysis (Chen et al., 2021b),
maximum likelihood estimation with hidden variables (Klusowski et al., 2019).
More importantly, learning by deep neural networks (DNN) is a wide class of
highly non-convex learning problems since the repeated composition of convex
functions is not convex.

We focus on the robust DNN regressions solved by stochastic gradient
algorithms (SGDs) in this paper, while the study in the future works on other
non-convex problems that need to use other algorithms (e.g., convex-concave-
procedure algorithms and expectation-maximization algorithms). By Theorem 2,
we can derive a general theorem (see Theorem 6 below) for DNN regressions.

We consider the DNN as the following compositional function class of functions
$f_{\theta}({x}):\mathbb{R}^{p}\to\mathbb{R}$

$\mathcal{N}\mathcal{N}(N,L):=\big{\\{}f_{\theta}({x})=W_{L}\sigma_{L}\left(W_{L-1}\sigma_{L-1}\left(\ldots
W_{1}\sigma_{1}\left(W_{0}{x}\right)\right)\right)|\,\theta:=(W_{0},\ldots,W_{L})\big{\\}},$
(16)

where $W_{l}\in\mathbb{R}^{N_{l+1}\times N_{l}}$ for $l=0,1,\dots,L$ with
$N_{0}=p$. Here $L$ represents the depth of this class of DNNs, each
activation function
$\sigma_{l}:\mathbb{R}^{N_{l}}\to\mathbb{R}^{N_{l}},l=0,1,2,\dots,L$, and
$\theta$ is the vectorized parameter consisting of weighted matrices with the
width $N=\max\left\\{N_{1},\ldots,N_{L}\right\\}$; see Fan et al. (2021) for
details.

For i.i.d. observations $\\{(X_{i},Y_{i})\\}_{i=1}^{n}$ and the given loss
function $l(\cdot,\cdot)$, we have risk function
$R_{l}(f_{\theta}):=\mathrm{E}[l(Y_{1},f_{\theta}({X_{1}}))]$ with a given
non-parametric function $f_{\theta}(\cdot)$. We assume that the parameter
$\theta^{*}$ of DNN is defined by

$\theta^{*}:=\mathop{\rm{argmin}}\limits_{\theta\in\Theta}R_{l}(f_{\theta}),~{}~{}f_{\theta}\in\mathcal{N}\mathcal{N}(N,L),$
(17)

where $\mathcal{N}\mathcal{N}(N,L)$ consists of feedforward neural networks
(16) parameterized by $\Theta$. The robust estimator of $\theta^{*}$ is
defined as

$\hat{\theta}_{n}:=\mathop{\arg\min}\limits_{\theta\in{\Theta}}\frac{1}{{n\alpha}}\sum_{i=1}^{n}\psi_{\lambda}\big{(}\alpha
l(Y_{i},f_{\theta}(X_{i}))\big{)}.$ (18)

Note that a typical way of estimating $\theta^{*}$ is by ERM estimator, i.e.
$\bar{\theta}_{n}={\rm argmin}_{\beta\in\Theta}\hat{R}^{n}_{l}(f_{\theta})$
with
$\hat{R}^{n}_{l}(f_{\theta})=\frac{1}{n}\sum_{i=1}^{n}l(Y_{i},f_{\theta}(X_{i}))$;
see Jiao et al. (2021); Shen et al. (2021b). However, ERM estimator is not
robust even when $f_{\theta}(\cdot)$ is linear.

Let
$\lambda(x)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon},~{}{\varepsilon}\in(0,1]$.
For a fixed $L$, Theorem 2 and Lipschitz property of DNNs [see Proposition 6
in Taheri et al. (2021)] imply

###### Theorem 6

Let ${\theta^{*}}$ be defined by (17) and $\hat{\theta}_{n}$ be given by (18)
under (C.1), where $f_{\theta}(\cdot)\in\mathcal{NN}(N,L)$ is the DNN function
indexed by parameter

$\theta\in\Theta:=\\{\theta:=(W_{1},\ldots,W_{L}):\|\theta\|_{\mathrm{F}}\leq
K~{}\text{and}~{}\max_{0\leq l\leq
L}{\sigma_{\max}}({W_{l}})\leq{M}\\}\subseteq\mathbb{R}^{\sum_{l=0}^{L}N_{l+1}N_{l}}.$

Suppose that $l(\cdot,\cdot)$ satisfies Lipschitz condition with a function
$D(x,y)$:

$|l(y,f_{\theta_{2}}(x))-l(y,f_{\theta_{1}}(x))|\leq
D(x,y){{|f_{\theta_{2}}(x)-f_{\theta_{1}}(x)|}},~{}\theta_{1},~{}\theta_{2}\in\Theta$.

If
$\alpha\leq\frac{1}{n^{1/({1+\varepsilon})}}(\frac{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}{[(2^{\varepsilon}+1)\sup_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]})^{1/({1+\varepsilon})}$ with
${\sup_{\theta\in\Theta}}R_{\lambda\circ l}(\theta)<\infty$, and
${{\rm{E}}}|\|X\|_{2}D(X,Y)|^{1+\varepsilon}<\infty$, then with probability at
least $1-2\delta$ one has

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{4(aM)^{L}\sqrt{L}}{n}{{\rm{E}}}[\|X\|_{2}D(X,Y)]+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2Kn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}$
$\displaystyle\times\left[\frac{2({2^{\varepsilon}}+1)}{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{-{({1+\varepsilon})^{-1}}}}+\frac{{{2^{\varepsilon}(2(aM)^{L}\sqrt{L})^{1+\varepsilon}}{{\rm{E}}}|\|X\|_{2}D(X,Y)|^{1+\varepsilon}}}{{(1+\varepsilon){n^{\varepsilon+1}}{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right].$ (19)

Based on (6), if the depth of DNNs is increasing with $n$, i.e. $L={L_{n}}$,
it allows
$\frac{{{{(aM)}^{L_{n}}}\sqrt{L_{n}}}}{n}\leq\frac{C}{{{n^{\varepsilon/(\varepsilon+1)}}}}$.
Then the _depth-sample condition_

${\rm{log}}(aM)\cdot{L_{n}}+\frac{1}{2}\log{L_{n}}\leq\log
C+\frac{1}{{\varepsilon+1}}\log n$

guarantees the rate $O({{{({{{p\log
n}}/{n}})}^{\varepsilon/(1+\varepsilon)}}})$ excess risk.

We finish this section by giving three concrete examples of robust DNN
regressions, and the models therein will be used in simulations or real data
studies.

######  Example 3 (Robust DNN LAD regression)

Suppose that i.i.d. observations
$\\{(Y_{i},X_{i})\sim(Y,X)\\}_{i=1}^{n}\in\mathbb{R}\times\mathbb{R}^{p}$
satisfies

$Y_{i}=f_{\theta^{*}}(X_{i})+e_{i},\quad\mathrm{E}(e_{i}|X_{i})=0,~{}i=1,2,\cdots,n,$

where $\\{e_{i}\sim e\\}_{i=1}^{n}$ is i.i.d. noise such that
$\mathrm{E}|e|^{1+\varepsilon}<\infty$ and
${\rm{E}}\left\|X\right\|_{2}^{1+\varepsilon}<\infty$. Similar to QR with
${\tau}=0.5$, the robust DNN LAD regression estimator (18) has loss function
$l(x,y,\theta)=|y-f_{\theta}(x)|$ and thus ${{R}_{\lambda\circ
l}}(\theta):=\frac{1}{1+\varepsilon}\mathrm{E}|Y-f_{\theta}(X)|^{1+\varepsilon}=\mathrm{E}|e|^{1+\varepsilon}/(1+\varepsilon)$
and we have $\sup_{\theta\in\Theta}{{R}_{\lambda\circ l}}(\theta)<\infty$ in
Theorem 6. Then (6) implies with probability at least $1-2\delta$

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{{4{{(aM)}^{L}}\sqrt{L}}}{n}{\rm{E}}{\left\|X\right\|_{2}}+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2Kn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}$
$\displaystyle\times\left[{{2({2^{\varepsilon}}+1){\left(\frac{\mathrm{E}|e|^{1+\varepsilon}}{1+\varepsilon}\right)^{{{(1+\varepsilon)}^{-1}}}}}}+\frac{{{2^{1+2\varepsilon}}{{(aM)}^{L(1+\varepsilon)}L^{(1+\varepsilon)/2}}{\rm{E}}\left\|X\right\|_{2}^{1+\varepsilon}}}{{{(\varepsilon+1)n^{\varepsilon+1}}{{[\mathrm{E}|e|^{1+\varepsilon}/(1+\varepsilon)]}^{\varepsilon/(\varepsilon+1)}}}}\right].$

######  Example 4 (Robust DNN logistic regression)

Assume that i.i.d. observations
$\\{(Y_{i},X_{i})\sim(Y,X)\\}_{i=1}^{n}\in\\{0,1\\}\times\mathbb{R}^{p}$
satisfy

$P({Y_{i}}=1|{X_{i}}):=\frac{{{e^{f_{\theta^{*}}(X_{i})}}}}{{1+{e^{f_{\theta^{*}}(X_{i})}}}},~{}~{}P({Y_{i}}=0|{X_{i}})=1-P({Y_{i}}=1|{X_{i}}),~{}f_{\theta^{*}}\in\mathcal{N}\mathcal{N}(N,L)$
(20)

with the moment conditions:
${\rm{E}}\left\|X\right\|_{2}^{1+\varepsilon}<\infty$ and
${\sup\limits_{\theta\in\Theta}}\mathrm{E}|{f_{\theta}(X)}|^{1+\varepsilon}<\infty$.
From (6), the robust DNN logistic regression estimator (18) has loss function
$l(x,y,\theta)=-[y{f_{\theta}(x)}-\log(1+{e^{f_{\theta}(x)}})]$. (6) shows
with probability at least $1-2\delta$

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{8(aM)^{L}\sqrt{L}}{n}{{\rm{E}}}\|X\|_{2}+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2Kn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}$
$\displaystyle\times\left[\frac{2({2^{\varepsilon}}+1)}{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{-{({1+\varepsilon})^{-1}}}}+\frac{{{2^{2+3\varepsilon}}{{(aM)}^{L(1+\varepsilon)}L^{(1+\varepsilon)/2}}{{\rm{E}}}\|X\|_{2}^{1+\varepsilon}}}{{(1+\varepsilon){n^{\varepsilon+1}}{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right],$

where ${\sup\limits_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta):=\frac{1}{1+\varepsilon}{\sup_{\theta\in\Theta}}\mathrm{E}[Y{f_{\theta}(X)}-\log(1+{e^{f_{\theta}(X)}})]^{1+\varepsilon}\leq
C_{1}{\sup\limits_{\theta\in\Theta}}\mathrm{E}|{f_{\theta}(X)}|^{1+\varepsilon}<\infty$
for constant $C_{1}>0$.

######  Example 5 (Robust DNN NBR)

Let $\eta>0$ be the known dispersion parameter. Suppose that i.i.d.
observations
$\\{(Y_{i},X_{i})\sim(Y,X)\\}_{i=1}^{n}\in\mathbb{Z}\times\mathbb{R}^{p}$
satisfy

$P({Y_{i}=y}|X_{i})=\frac{\Gamma(\eta+y)}{\Gamma(\eta)y!}(\frac{\mu_{i}}{\eta+\mu_{i}})^{y}(\frac{\eta}{\eta+\mu_{i}})^{\eta},~{}\text{with}~{}\log\mu_{i}=f_{\theta^{*}}(X_{i}),~{}f_{\theta^{*}}\in\mathcal{N}\mathcal{N}(N,L)$

with the moment conditions:
${\rm{E}}[\|X\|_{2}(Y+\eta)]^{1+\varepsilon}<\infty$ and
${\sup\limits_{\theta\in\Theta}}{{\rm{E}}}|Y{f_{\theta}(X)}|^{1+\varepsilon}<\infty$.
From (6), the robust DNN NBR estimator (18) has loss function
$l(x,y,\theta)=-y[{f_{\theta}(x)}-\log(\eta+{e^{f_{\theta}(x)}})]-\eta\log(\eta+{e^{f_{\theta}(x)}})$.
The (6) shows with probability at least $1-2\delta$,

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\beta^{*}})$
$\displaystyle\leq\frac{4(aM)^{L}\sqrt{L}}{n}{{\rm{E}}}[\|X\|_{2}(Y+\eta)]+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2Kn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}$
$\displaystyle\times\left[\frac{2({2^{\varepsilon}}+1)}{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{-{({1+\varepsilon})^{-1}}}}+\frac{{{2^{1+2\varepsilon}}{{(aM)}^{L(1+\varepsilon)}L^{(1+\varepsilon)/2}}{{\rm{E}}}[\|X\|_{2}(Y+\eta)]^{1+\varepsilon}}}{{(1+\varepsilon){n^{\varepsilon+1}}{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right],$

where ${\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta):=\frac{1}{1+\varepsilon}{\sup_{\theta\in\Theta}}{{\rm{E}}}{\\{Y[{f_{\theta}(X)}-\log(\eta+{e^{f_{\theta}(X)}})]-\eta\log(\eta+{e^{f_{\theta}(X)}})\\}^{1+\varepsilon}}\leq
C_{1}{\sup\limits_{\theta\in\Theta}}{{\rm{E}}}|Y{f_{\theta}(X)}|^{1+\varepsilon}+C_{2}{\sup\limits_{\theta\in\Theta}}{{\rm{E}}}|{f_{\theta}(X)}|^{1+\varepsilon}<\infty$
for constants $C_{1},C_{2}>0$.

##  4 Simulation and real data studies

###  4.1 Simulation

####  4.1.1 Normal regression models

We illustrate the effectiveness of regressions based on log-truncated ERM by
the numerical experiments of standard logistic regression and negative
binomial regression. Due to the condition (C.0): $\|\theta\|_{2}\leq r$, the
regression problem (5) is equivalent to the following $\ell_{2}$-regularized
optimization:

$\hat{\theta}_{n}(\alpha,\rho):=\mathop{\arg\min}\limits_{\theta\in{\Theta}}\\{{{\hat{R}}_{\psi,l,\alpha}}(\theta)+\rho\|\theta\|_{2}\\},$
(21)

where
${{\hat{R}}_{\psi,l,\alpha}}(\theta):=\frac{1}{{n\alpha}}\sum_{i=1}^{n}\psi_{\lambda}[\alpha
l({Y_{i}},{{X}}_{i},\theta)]$ with $l(y,x,\theta)$ being logistic loss or
negative binomial loss, $\rho>0$ is the penalty parameter for
$\ell_{2}$-regularization, and $\alpha>0$ is another tuning parameter to be
chosen. In practice, this optimization problem is solved by stochastic
gradient descent (SGD) as the following:

$\theta_{t+1}=\theta_{t}-\frac{r_{t}}{\alpha}\nabla_{\theta}\\{\psi_{\lambda}[\alpha
l(Y_{i_{t}},X_{i_{t}},\theta_{t})]+\rho\|\theta_{t}\|_{2}\\},t=0,1,2,\cdots,$
(22)

where $i_{t}$ denotes a random sampled index, $\\{r_{t}\\}$ is the learning
rate. We take binary search to find the optimal parameter pair $(\alpha,\rho)$
in certain effective subset of $\mathbb{R}^{+}\times\mathbb{R}^{+}$.

For the regressions without truncation, the corresponding optimization problem
is

$\hat{\theta}_{n}(\rho):=\mathop{\arg\min}\limits_{\theta\in{\Theta}}\\{\frac{1}{{n}}\sum_{i=1}^{n}l({Y_{i}},{{X}}_{i},\theta)+\rho\|\theta\|_{2}\\},$
(23)

where $\rho$ is the penalty parameter for regularization, this optimization
problem can be solved by the following SGD:

$\theta_{t+1}=\theta_{t}-r_{t}\nabla_{\theta}\\{l(Y_{i_{t}},X_{i_{t}},\theta_{t})+\rho\|\theta_{t}\|_{2}\\},t=0,1,2,\cdots.$
(24)

Here we also take binary search to find the optimal parameter $\rho$ in
certain effective subset of $\mathbb{R}^{+}$.

For $\mathbb{R}^{p}$ valued covariates $X_{1},...,X_{n}$, for each $X_{i}$, it
can be written as

$X_{i}=X_{i}^{\prime}+\xi_{i},$

where $X_{1}^{\prime},...,X_{n}^{\prime}$ are i.i.d. $\mathbb{R}^{p}$-valued
random vectors with normal distribution $N({\bf 0},{\bf I}_{p\times p})$, and
$\xi_{1},...,\xi_{n}$ are i.i.d. $\mathbb{R}^{p}$-valued random vectors which
satisfy one of the following conditions:

  * (1)

Pareto noise (heavy tail): the noise $\xi_{i}=(\xi_{i1},...,\xi_{ip})^{\top}$
whose entries $\xi_{ij}$ are independently drawn from Pareto distribution with
shape parameter $\beta$, we will take $\beta\in\\{1.6,1.8,2.01,4.01,6.01\\}$.

  * (2)

Uniform noise (outlier): the noise $\xi_{i}=Z\xi_{i}^{\prime}$ where
$\xi_{i}^{\prime}=(\xi_{i1}^{\prime},...,\xi_{ip}^{\prime})^{\top}$ are
independently drawn from uniform distribution $U(10,20)$ in logistic
regression and $U(1,5)$ in negative binomial regression (with dispersion
parameter $\eta=0.1$), and $Z={\rm{diag}}(\zeta_{1},...,\zeta_{p})$ where
$\zeta_{1},...,\zeta_{p}$ are i.i.d. Bernoulli distribution with
$\psi\in(0,1)$. We will choose $\psi\in\\{0.3,0.5,0.8\\}$.

For each fixed pair $(p,n)$ in simulation, the true value
$\theta=(\theta_{1},...,\theta_{p})^{\top}\in\mathbb{R}^{p}$ is designed by
drawing each $\theta_{j}$ from $U(0,1)$ independently. We choose log-truncated
function as (4) with $\lambda(|x|)=\frac{|x|^{1+\varepsilon}}{1+\varepsilon}$.
When $\\{\xi_{i}\\}$ are Pareto noises, we choose $1+\varepsilon=1.5$ as
$\beta\in\\{1.6,1.8\\}$ and $1+\varepsilon=2.0$ as
$\beta\in\\{2.01,4.01,6.01\\}$. When $\xi_{i}$ are uniform noises, we always
choose $1+\varepsilon=2.0$.

We will conduct experiments for the following three cases:

$(p,n)=(100,200),(200,500),(1000,1000).$

Table $\ref{tab:LGPareto}-$3 present the comparison of average
$\ell_{2}$-estimation errors and standard errors (in bracket) for predicted
coefficients of logistic regression with 100 replications. The log-truncated
estimators perform much better than standard regression estimators when Pareto
noises contaminate covariates; the larger the sample size $n$ is, the more
easily the effectiveness of the log-truncated estimators can be seen. We can
see the same advantages of the log-truncated estimators from the table below
for the uniform noise model. The results of simulations for NBR are displayed
in the appendix.

####  4.1.2 Deep regression models

We also apply the deep logistic regression and deep negative binomial
regression in Example $2$ and $3$ to the simulated data. We consider the non-
parametric model $f_{\theta^{\ast}}(X_{i})$ for $i=1,...,n$ in (20) and (5)
satisfying

$f_{\theta^{\ast}}(X_{i})=f_{0}(X_{i}+\xi_{i}),\ i=1,...,n,$ (25)

where $\\{X_{i}\\}_{i}^{n}$ are i.i.d. $\mathbb{R}^{p}$-valued random vectors
with normal distribution $N({\bf 0},{\bf I}_{p\times p})$,
$\\{\xi_{i}\\}_{i}^{n}\in\mathbb{R}^{p}$ also follow the above noise
condition, Pareto and uniform noise. For the Pareto noise, the shape parameter
$\beta$ are set $\\{1.6,1.8,2.01,2.80,3.01\\}$. For the uniform noise, we
choose the probability $\psi$ in Bernoulli distribution in
$\\{0.3,0.5,0.8\\}$.

The continuous function $f_{0}(X_{i})$: $\mathbb{R}^{p}\to\mathbb{R}$,
subjects to the additive model for each $X_{i},i=1,...,n$:

$f_{0}(X_{i})=\sum_{j=1}^{p}f_{j}(X_{ij}),i=1,...,n,$

where $f_{j}$ are smooth functions selected in a set of univariate functions

$\mathcal{G}=\\{g_{k}(x),k=1,2,3\\},$

with $g_{1}(x)=-10x+0.3$, $g_{2}(x)=0.7x^{3}-0.2x^{2}+0.3x-0.3$,
$g_{3}(x)=0.3\sin(x)\sqrt{|x|}$ referring to the experiments in Shen et al.
(2021b). For the choice of $f_{j},j=1,...,p$, we first uniformly generate a
random matrix with $n\times p$ integers in $\\{1,2,3\\}$. Then, viewing the
integers in generated random matrix as the indices, we choose $f_{j}$ in
$\mathcal{G}$ for $j=1,...,p$.

For the network configuration, we minimize the ERM problem of logistic
regression by ReLU activated 2-layers DNN model and network width being
$(p,0.6p,0.4p,1)$. The batch size is set to be $n/4$ in each case. For DNN
logistic regression, we compute the average accuracy of predictors
$\\{\hat{Y}_{i}\\}_{i=1}^{n}$ of responses $\\{Y_{i}\\}_{i=1}^{n}\in\\{0,1\\}$
as $p=50,n=1000$, satisfying

$\text{Accuracy}=\frac{1}{n}\sum_{i=1}^{n}\mathbf{1}(\hat{Y}_{i}=Y_{i})\times
100\%.$

Table 4 shows the average accuracy for DNN logistic regression under Pareto
and uniform noise conditions with 100 replications. Different values of
$1+\varepsilon$ are selected, see Table 4. The results reveal that the log-
truncated DNN logistic regression has higher flexibility and stronger
robustness than standard DNN logistic regression for fitting the heavy-tail
data and the data with outliers. The results of DNN NBR are demonstrated in
the appendix.

Table 2: Comparison of average $\ell_{2}$-estimation error for logistic
regression on Pareto noise model.

$\ell_{2}$-estimation error for logistic regression
---
$p$ = 100, $n$ = 200
$1+\varepsilon$ |  | $\beta$ |  | Truncation |  | Non-truncation
1.5 |  | 1.60 |  | 0.1152(0.0141) |  | 0.1270(0.0170)
1.5 |  | 1.80 |  | 0.0971(0.0997) |  | 0.1167(0.0130)
2.0 |  | 2.01 |  | 0.1189(0.0164) |  | 0.1931(0.0200)
2.0 |  | 4.01 |  | 0.1128(0.0153) |  | 0.1531(0.0157)
2.0 |  | 6.01 |  | 0.1078(0.0131) |  | 0.1411(0.0158)
$p$ = 200, $n$ = 500
$1+\varepsilon$ |  | $\beta$ |  | Truncation |  | Non-truncation
1.5 |  | 1.60 |  | 0.1070(0.0105) |  | 0.1179(0.0102)
1.5 |  | 1.80 |  | 0.0944(0.0092) |  | 0.1183(0.0010)
2.0 |  | 2.01 |  | 0.1065(0.0099) |  | 0.1870(0.0080)
2.0 |  | 4.01 |  | 0.0959(0.0085) |  | 0.1485(0.0078)
2.0 |  | 6.01 |  | 0.0931(0.0106) |  | 0.1388(0.0074)
$p$ = 1000, $n$ = 1000
$1+\varepsilon$ |  | $\beta$ |  | Truncation |  | Non-truncation
1.5 |  | 1.60 |  | 0.1233(0.0052) |  | 0.1471(0.0052)
1.5 |  | 1.80 |  | 0.1222(0.0044) |  | 0.1462(0.0055)
2.0 |  | 2.01 |  | 0.1345(0.0053) |  | 0.2495(0.0041)
2.0 |  | 4.01 |  | 0.1301(0.0046) |  | 0.2372(0.0040)
2.0 |  | 6.01 |  | 0.1302(0.0042) |  | 0.2356(0.0039)

Table 3: Comparison of average $\ell_{2}$-estimation error for logistic
regression on Uniform noise model.

$\ell_{2}$-estimation error for logistic regression
---
$p$ = 100, $n$ = 200
$1+\varepsilon$ |  | $\psi$ |  | Truncation |  | Non-truncation
2.0 |  | 0.3 |  | 0.1189(0.0164) |  | 0.3156(0.0148)
2.0 |  | 0.5 |  | 0.1283(0.0176) |  | 0.2947(0.0209)
2.0 |  | 0.8 |  | 0.1517(0.0182) |  | 0.2631(0.0228)
$p$ = 200, $n$ = 500
$1+\varepsilon$ |  | $\psi$ |  | Truncation |  | Non-truncation
2.0 |  | 0.3 |  | 0.1114(0.0105) |  | 0.3106(0.0100)
2.0 |  | 0.5 |  | 0.1214(0.0111) |  | 0.2837(0.0118)
2.0 |  | 0.8 |  | 0.1405(0.0117) |  | 0.2489(0.0113)
$p$ = 1000, $n$ = 1000
$1+\varepsilon$ |  | $\psi$ |  | Truncation |  | Non-truncation
2.0 |  | 0.3 |  | 0.1397(0.0051) |  | 0.3322(0.0034)
2.0 |  | 0.5 |  | 0.1369(0.0049) |  | 0.3129(0.0040)
2.0 |  | 0.8 |  | 0.1413(0.0050) |  | 0.2904(0.0039)

Table 4: Comparison of average accuracy for DNN logistic regression under two
noise conditions.

Accuracy (%) for DNN logistic regression
---
$p$ = 500, $n$ = 1000
$1+\varepsilon$ |  | Pareto ($\beta$) |  | Truncation |  | Non-truncation
1.3 |  | 1.80 |  | 93.12(0.0269) |  | 92.66(0.0334)
1.7 |  | 1.90 |  | 86.77(0.0009) |  | 82.94(0.0003)
2.0 |  | 2.01 |  | 81.91(0.0240) |  | 76.81(0.0455)
2.0 |  | 2.80 |  | 78.86(0.0010) |  | 74.89(0.0030)
2.0 |  | 3.01 |  | 85.19(0.0189) |  | 81.95(0.0451)
$1+\varepsilon$ |  | Uniform ($\psi$) |  | Truncation |  | Non-truncation
2.0 |  | 0.3 |  | 78.83(0.0041) |  | 60.57(0.0046)
2.0 |  | 0.5 |  | 73.31(0.0001) |  | 59.06(0.0215)
2.0 |  | 0.8 |  | 70.37(0.0007) |  | 53.05(0.0259)

###  4.2 Real data analysis

####  4.2.1 Boston housing dataset

We use the Boston housing dataset provided by the python library Scikit-Learn
to learn a log-truncated deep LAD model. Boston housing dataset contains
$n=506$ cases, and each case includes 14 variables. We aim to predict the
variable, Median Value (MEDV) of Owner-Occupied Housing Units as output, by
the remaining 13 ones. To this end, we randomly split the dataset into two
groups, one as the training set and the other as the testing set, and train a
$l_{2}$-regularized standard LAD regression in (21) and a 3-layers DNN LAD
regression model separately.

We denote the variable MEDV by $Y\in\mathbb{R}$ and the other 13 variables by
$X\in\mathbb{R}^{13}$, so the Boston housing dataset can be represented as
$\\{(Y_{i},X_{i}):i=1,2,...,506\\}$. In our experiment, the training set has
$n_{1}=339$ cases and the testing set has $n_{2}=167$ cases, denoted by
$\\{(Y_{tr,1},X_{tr,1}),...,(Y_{tr,n_{1}},X_{tr,n_{1}})\\}$ and
$\\{(Y_{te,1},X_{te,1}),...,(Y_{te,n_{2}},X_{te,n_{2}})\\}$ respectively. We
use the training set to train a 3-layers DNN model, and feed the testing set
$\\{X_{te,1},...,X_{te,n_{2}}\\}$ into the model to get a prediction set
$\\{\hat{Y}_{1},...,\hat{Y}_{n_{2}}\\}$. To assess the obtained DNN model, we
compute the absolute average errors (MAEs) of the prediction, which are
defined as

$\text{MAE}:=\frac{1}{n_{2}}\sum_{i=1}^{n_{2}}|\hat{Y}_{i}-Y_{te,i}|.$

We also use the above approach to train a $l_{2}$-regularized standard LAD
regression and compute the MAEs of the prediction separately. The high order
function in (4) is chosen as
$\lambda(|x|)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon}$. To select an
appropriate $1+\varepsilon$, we consider ten values of $1+\varepsilon$ in
$(1,2]$ with $\varepsilon=0.1,0.2,...,1$; see Table 5. Table 5 indicates that
the truncated $l_{2}$-regularized standard LAD regression has smaller
prediction errors than its non-truncated version as
$1+\varepsilon=1.1,1.2,1.8$. The truncated deep LAD model outperforms than
non-truncated deep LAD model on the whole. In particular, as
$1+\varepsilon=1.8$, i.e., $\varepsilon=0.8$, the prediction errors are
smallest for both truncated $l_{2}$-regularized standard LAD model and
truncated deep LAD model. Thus, we choose $1+\varepsilon=1.8$ for these two
models.

Table 5: Comparison of MAEs on Boston housing dataset.

| Standard LAD regression |  | Deep LAD regression
---|---|---|---
$1+\varepsilon$ | Truncation | Non-truncation |  | Truncation | Non-truncation
1.1 | 6.5203 | 6.5972 |  | 6.5115 | 6.9845
1.2 | 6.5464 | 6.5972 |  | 6.6729 | 6.9845
1.3 | 6.5517 | 6.5972 |  | 6.2602 | 6.9845
1.4 | 6.5737 | 6.5972 |  | 6.6777 | 6.9845
1.5 | 6.5756 | 6.5972 |  | 6.5351 | 6.9845
1.6 | 6.6196 | 6.5972 |  | 6.7033 | 6.9845
1.7 | 6.5853 | 6.5972 |  | 6.6189 | 6.9845
1.8 | 6.4041 | 6.5972 |  | 6.0884 | 6.9845
1.9 | 6.5727 | 6.5972 |  | 6.6772 | 6.9845
2.0 | 6.5492 | 6.5972 |  | 6.3412 | 6.9845

The optimal $\varepsilon$ in Table 5 can be roughly interpreted by Theorem 2,
from which we can see that there is a trade-off between the constant
${{\rm{E}}}|{H}(Y,X)|^{1+\varepsilon}$ and the excess risk convergence rate
$\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}$.

####  4.2.2 MNIST database

We use a handwritten digits database MNIST to learn a robust 6-layers DNN LAD
model compared with the standard 6-layers DNN LAD model, which is learned in
the same way. The activation function in (16) is ReLU. MNIST database contains
70000 28$\times$28 grayscale images of the 10 digits. In our experiment, we
treat the digit images as the input variables $X$ and their corresponding
labels as output $Y\in\\{0,1,2,\cdots,9\\}$. There are 60000 training images
and 10000 testing images in this database, denoted by
$\\{(X_{tr,i},Y_{tr,i})\\}\\}_{i=1}^{60000}$ and
$\\{(X_{te,i},Y_{te,i})\\}_{i=1}^{10000}$ respectively. Firstly, we normalize
the 28$\times$28 pixels of each image in range $(-1,1)$. Then, we set the
batch size to 64 and train the DNN over 60000 training images. To evaluate the
performance of the trained DNN model, we compute their classification accuracy
over 10000 testing images. That is:

$\text{Accuracy}:=\frac{1}{10000}\sum_{i=1}^{10000}\mathbf{1}(\hat{Y}_{te,i}=Y_{te,i})\times
100\%,$ (26)

where $\mathbf{1}(\cdot)$ denotes the indicator function,
$\\{\hat{Y}_{te,i}\\}_{i=1}^{10000}$ are the predictors from 10000 testing
images.

Three types of noises are added to the original data: Gaussian noise
$N(0.5,2)$, uniform noise $U(2,5)$ and Pareto noise with shape parameter
$\beta=2.01$ on the 28$\times$28 pixels. For each batch, we randomly select
$20\%,50\%$ and $80\%$ noises to contaminate the data. The truncation function
is also chosen as (4) with
$\lambda(|x|)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon}$. We selected the
optimal values of $1+\varepsilon$ in $(1,2]$ with the interval of 0.1 for the
three classes of noise. We repeat the experiment 100 times and record the
average accuracy, see Table 6 (the values in bracket are standard errors of
accuracy). The results in Table 6 show that the classification results of the
truncated DNN LAD model are better than the standard DNN LAD model. Especially
under strong noise disturbance, the truncated DNN LAD model is more robust
than the standard DNN LAD model.

We can see from the tables that the optimal index $1+\varepsilon$ does not
change according to the proportions of noises, and this is a significant
advantage of our truncation regression models.

Table 6: Comparison of classification accuracy on MNIST dataset Accuracy (%)
---
Gaussian Noise $N(0.5,2)$
$1+\varepsilon$ | Noise proportion | Truncation |  | Non-truncation
2.0 | 20% | 91.58 (0.17) |  | 90.62 (0.16)
2.0 | 50% | 84.23 (0.28) |  | 82.33 (0.29)
2.0 | 80% | 77.60 (0.33) |  | 74.50 (0.28)
Uniform Noise $U(2,5)$
$1+\varepsilon$ | Noise proportion | Truncation |  | Non-truncation
1.5 | 20% | 96.18 (0.09) |  | 95.75 (0.09)
1.5 | 50% | 94.76 (0.15) |  | 92.95 (0.12)
1.5 | 80% | 92.08 (0.19) |  | 84.04 (0.19)
Pareto Noise $\beta=2.01$
$1+\varepsilon$ | Noise proportion | Truncation |  | Non-truncation
1.3 | 20% | 88.96 (0.16) |  | 87.85 (0.20)
1.3 | 50% | 77.35 (0.33) |  | 72.99 (0.29)
1.3 | 80% | 66.32 (0.44) |  | 61.50 (0.34)

####  4.2.3 German health care demand database

To study count data regressions, we apply the robust DNN NBR to the German
health care demand database, provided by the German Socioeconomic Panel
(GSOEP), see Riphahn et al. (2003). The data source is

http://qed.econ.queensu.ca/jae/2003-v18.4/riphahn-wambach-million/.

This dataset contains 27326 patients (with identification number) observed
during 7 years $\\{1984,1985,1986,1987,1988,1991,1994\\}$, where sample sizes
for each year are

$\\{3874,3794,3792,3666,4483,4340,3377\\}$,

respectively. For a specific patient, there are 21 covariates and 2
alternatives, independent variables DOCVIS (number of doctor visits within the
last three months prior to the survey) and HOSPVIS (number of hospital visits
in the last calendar year) can be chosen. Table 7 gives the results of
descriptive statistics of the variables in the dataset. We notice that the
distributions of some covariates are heavy-tailed by showing a high
coefficient of kurtosis such as _age, handper, hhninc, and educ_. To this end,
our proposed robust DNN NBR is used to explore the relationship between
independent variables DOCVIS, HOSPVIS, and 21 covariates separately.

The problem is that the dispersion parameter $\eta$ in NBR is unknown. We need
to first estimate it in DNN NBR by maximizing the joint log-likelihood
function.

Next, the obtained $\hat{\eta}$ is plugged into the DNN NBR optimization in
Example 5. We randomly split the German health care demand dataset into
training dataset $(X_{train},Y_{train})$ including 33% observations and
testing dataset $(X_{test},Y_{test})$ including 67% observations for each
year. The truncation function is
$\lambda(|x|)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon}$. We consider
$1+\varepsilon$ in $(1,2]$ with steps 0.1 and select the optimal value of
$1+\varepsilon$. Firstly, we use the training dataset to train a 3-layers DNN,
and the activation function is softmax. Then, we compute the absolute average
errors (MAEs) of the predictive value of $Y_{test}$ as the measure of
accuracy. A standard DNN NBR is also trained in the same way. The MAEs of two
DNN NBRs are shown in Table 8. The results indicate that the truncated DNN NBR
produces smaller prediction errors than the standard DNN NBR under the partly
heavy-tailed inputs.

Table 7: Descriptive statistics of German health care demand database

Variables | Description | Mean | SD | Skew | Kurtosis
---|---|---|---|---|---
docvis | number of doctor visits in last three months | 3.183525 | 5.68969 | 5.184572 | 47.09992
hospvis | number of hospital visits in last calendar year | 0.138257 | 0.884339 | 26.66923 | 1098.842
age | age in years | 43.52569 | 11.33025 | 0.084424 | -1.15198
hsat | health satisfaction, coded 0 (low) - 10 (high) | 6.785426 | 2.293715 | -0.69139 | 0.128572
handdum | handicapped = 1; otherwise = 0 | 0.214046 | 0.410056 | 1.394705 | -0.05433
handper | degree of handicap in percent (0 - 100) | 7.012285 | 19.26463 | 3.149664 | 9.345931
hhninc | household nominal monthly net income in German marks / 1000 | 3.520836 | 1.769083 | 2.426049 | 17.79162
hhkids | children under age 16 in the household = 1; otherwise = 0 | 0.40273 | 0.490456 | 0.396637 | -1.84275
educ | years of schooling | 11.32063 | 2.324885 | 1.622374 | 2.286639
married | married = 1; otherwise = 0 | 0.758618 | 0.427929 | -1.20865 | -0.53918
haupts | highest schooling degree is Hauptschul degree = 1; otherwise = 0 | 0.624277 | 0.484318 | -0.51319 | -1.7367
reals | highest schooling degree is Realschul degree = 1; otherwise = 0 | 0.196809 | 0.397594 | 1.525074 | 0.325861
fachhs | highest schooling degree is Polytechnical degree = 1; otherwise = 0 | 0.04084 | 0.197924 | 4.639599 | 19.52659
abitur | highest schooling degree is Abitur = 1; otherwise = 0 | 0.117031 | 0.321464 | 2.382569 | 3.67677
univ | highest schooling degree is university degree = 1; otherwise = 0 | 0.071946 | 0.258403 | 3.312945 | 8.975935
working | employed = 1; otherwise = 0 | 0.677048 | 0.467613 | -0.75721 | -1.42668
bluec | blue collar employee = 1; otherwise = 0 | 0.243761 | 0.429358 | 1.193551 | -0.57546
whitec | white collar employee = 1; otherwise = 0 | 0.299605 | 0.458093 | 0.874878 | -1.23463
self | self employed = 1; otherwise = 0 | 0.062175 | 0.241478 | 3.626077 | 11.14884
beamt | civil servant = 1; otherwise = 0 | 0.074691 | 0.262897 | 3.235445 | 8.468415
public | insured in public health insurance = 1; otherwise = 0 | 0.885713 | 0.318165 | -2.42452 | 3.87845
addon | insured by add-on insurance = 1; otherswise = 0 | 0.01881 | 0.135856 | 7.083579 | 48.17885
year | calendar year of the observation (1984, 1985, 1986, 1987, 1988, 1991, 1994) |  |  |  |

Table 8: Comparison of MAEs on German health care demand dataset MAEs for DNN
NBR
---
|  | DOCVIS |  | HOSPVIS
Year | $1+\varepsilon$ | Truncation |  | Non-truncation |  | Truncation |  | Non-truncation
1984 | 2.0 | 3.5265 |  | 4.1360 |  | 0.4895 |  | 0.6132
1985 | 1.5 | 4.3428 |  | 5.5667 |  | 2.3165 |  | 2.3441
1986 | 2.0 | 4.3369 |  | 4.2114 |  | 0.6446 |  | 0.8558
1987 | 1.5 | 3.9456 |  | 4.0242 |  | 1.7162 |  | 1.8948
1988 | 1.5 | 2.5847 |  | 2.5697 |  | 0.4676 |  | 0.6100
1991 | 1.5 | 2.4537 |  | 2.7494 |  | 0.8253 |  | 0.9670
1994 | 1.5 | 3.4578 |  | 3.6456 |  | 0.6861 |  | 0.6971

##  5 Discussion and Future Study

From the simulations and real data analysis above, we can see that the
selection of the $1+\varepsilon$ is a very important issue for prediction. The
tail index estimation has been intensively studied in extreme-value
statistics; see Fedotenkov (2020) for a review. We leave the research of
estimating the index $1+\varepsilon$ to the future study.

As mentioned above, Theorem 2 can be applied to study other robust statistical
learning models, in which one may have to design specific algorithms according
to concrete problems at hand. We conclude this paper by the following two
examples, robust two-component mixed linear regression and robust non-negative
matrix factorization, whose algorithms are different from SGD. Here we only
roughly address their theoretical results and leave detailed study to the
future work.

Robust two-component mixed linear regression. One challenging non-convex
problem is the mixture of two linear regressions. Suppose
$\\{(Y_{i},X_{i})\sim(Y,X)\\}_{i=1}^{n}$ are
$\mathbb{R}\times\mathbb{R}^{p}$-valued i.i.d. random variables. With
probability $\pi$, $(X,Y)$ has conditional density function
$p(y,{x^{\top}}{\eta_{0}}),~{}{\eta_{0}}\in{\mathbb{R}^{p}}$ for $Y=y|X=x$,
and with probability $1-\pi$, $(X,Y)$ has conditional density function
$p(y,{x^{\top}}{\eta_{1}}),~{}{\eta_{1}}\in{\mathbb{R}^{p}}$ for $Y=y|X=x$,
where ${\eta_{0}}$ and ${\eta_{1}}$ are unknown coefficients. Given the input
$x$, the negative log-likelihood function of the output $y$ is

$\displaystyle l(y,x,;{\pi,\eta_{0},\eta_{1}})=-\log[\pi
p(y,{x^{\top}}{\eta_{0}})+(1-\pi)p(y,{x^{\top}}{\eta_{1}})],$ (27)

where $\pi\in(0,1)$ is an unknown mixing probability. Mei et al. (2018)
studied Gaussian mixture models, while Khamaru and Wainwright (2019)
considered mixture density estimation under subexponential condition.

In order to fit this example to our theory, we write
$\theta:=\left({\pi,{\eta_{0}},{\eta_{1}}}\right)$. Under the following moment
conditions: ${\rm{E}}({\left\|X\right\|_{2}}|Y|)^{1+\varepsilon}<\infty$,
${\rm{E}}\left\|X\right\|_{2}^{2+2\varepsilon}<\infty$, by Theorem 2, we can
obtain an excess risk in the order of $O({{{({{{p\log
n}}/{n}})}^{\varepsilon/(1+\varepsilon)}}})$.

Robust non-negative matrix factorization. Given $n$ vector samples
$\\{x_{i}\\}_{i=1}^{n}\in\mathbb{R}^{p}$ arranged in a nonnegative matrix
$\mathrm{S}:=\left[x_{1},\ldots,x_{n}\right]\in\mathbb{R}_{+}^{p\times n}$ and
positive integer $r\geq 1$, the log-truncated non-negative matrix
factorization (NMF) considers the decomposition for $\mathrm{S}$ into a
product of two non-negative matrices:

$\mathrm{S}=\mathrm{B}\mathrm{C}+\mathrm{R},$

where $\mathrm{B}\in\mathbb{R}_{+}^{p\times r}$ is the given basis,
$\mathrm{C}=\left[h_{1},\ldots,h_{n}\right]\in\mathbb{R}_{+}^{r\times n}$ is
the coefficients, and $\mathrm{R}\in\mathbb{R}^{p\times n}$ is the random
error matrix.

By minimizing the $\ell_{2}$-distance between their product and the original
data matrix, the ordinary NMF decomposes a data matrix into the product of two
lower dimensional non-negative factor matrices $\mathrm{B}$ and $\mathrm{C}$.
When the original data matrix is corrupted by heavy-tailed outliers that
seriously violate the second moment assumption (Guan et al., 2017). We assume
the element of $R_{ij}$ in $\mathrm{R}$ has only $(1+\varepsilon)$-th moment,
i.e. $\sup_{i\in[n],j\in[p]}{\rm{E}}|R_{ij}|^{1+\varepsilon}<\infty$, it is of
interest to study the log-truncated NMF and the robust algorithm by using the
following non-convex optimization problem:

$\hat{\mathrm{C}}=\mathop{\arg\min}\limits_{\mathrm{C}\in\mathbb{R}_{+}^{r\times
n}}\frac{1}{{np\alpha}}\sum\limits_{i\in[n],j\in[p]}{\psi_{\lambda}}\left[\alpha\left({\mathrm{S}-\mathrm{B}\mathrm{C}}\right)_{ij}^{2}\right]$

where
$\lambda(x)=\frac{1}{1+\varepsilon}|x|^{1+\varepsilon},~{}{\varepsilon}\in(0,1]$
and the $\alpha$ is tuning parameter implicitly determined by the random error
matrix. We expect that $\hat{\mathrm{C}}$ is able to learn a subspace on a
dataset though the original data matrix that is contaminated by a heavy-tailed
noise matrix.

Acknowledgments

L. Xu is supported in part by NSFC grant (No.12071499), Macao S.A.R grant FDCT
0090/2019/A2 and University of Macau grant MYRG2018-00133-FST. F. Yao is
supported by National Natural Science Foundation of China Grants No.11931001
and 11871080, the LMAM, and the Key Laboratory of Mathematical Economics and
Quantitative Finance (Peking University), Ministry of Education, National Key
R&D Program of China Grant (No. 2020YFE0204200). H. Zhang is supported in part
by the University of Macau under UM Macao Talent Programme (UMMTP-2020-01) and
National Natural Science Foundation of China Grants (No.12101630).

## References

  * Baraud et al. (2017) Yannick Baraud, Lucien Birgé, and Mathieu Sart.  A new method for estimation and model selection: $rho$-estimation.  _Inventiones mathematicae_ , 207(2):425–517, 2017\.
  * Bartlett and Mendelson (2006) Peter L Bartlett and Shahar Mendelson.  Empirical minimization.  _Probability theory and related fields_ , 135(3):311–334, 2006.
  * Brownlees et al. (2015) Christian Brownlees, Emilien Joly, Gábor Lugosi, et al.  Empirical risk minimization for heavy-tailed losses.  _Annals of Statistics_ , 43(6):2507–2536, 2015\.
  * Catoni (2012) Olivier Catoni.  Challenging the empirical mean and empirical variance: a deviation study.  In _Annales de l’IHP Probabilités et statistiques_ , volume 48, pages 1148–1185, 2012.
  * Chen et al. (2021a) Peng Chen, Xinghu Jin, Xiang Li, and Lihu Xu.  A generalized catoni’s m-estimator under finite $\alpha$-th moment assumption with $\alpha\in(1,2)$.  _Electronic Journal of Statistics_ , 15(2):5523–5544, 2021a.
  * Chen et al. (2021b) Yuxin Chen, Jianqing Fan, Cong Ma, and Yuling Yan.  Bridging convex and nonconvex optimization in robust pca: Noise, outliers and missing data.  _The Annals of Statistics_ , 49(5):2948–2971, 2021b.
  * Chi (2010) Zhiyi Chi.  A local stochastic lipschitz condition with application to lasso for high dimensional generalized linear models.  _arXiv preprint arXiv:1009.1052_ , 2010.
  * Chinot et al. (2019) Geoffrey Chinot, Guillaume Lecué, and Matthieu Lerasle.  Robust statistical learning with lipschitz and convex loss functions.  _Probability Theory and related fields_ , pages 1–44, 2019.
  * Fan et al. (2017) Jianqing Fan, Quefeng Li, and Yuyan Wang.  Estimation of high dimensional mean regression in the absence of symmetry and light tail assumptions.  _Journal of the Royal Statistical Society. Series B, Statistical methodology_ , 79(1):247, 2017.
  * Fan et al. (2021) Jianqing Fan, Cong Ma, and Yiqiao Zhong.  A selective overview of deep learning.  _Statistical Science_ , 36(2):264–290, 2021.
  * Fedotenkov (2020) Igor Fedotenkov.  A review of more than one hundred pareto-tail index estimators.  _Statistica_ , 80(3):245–299, 2020.
  * Guan et al. (2017) Naiyang Guan, Tongliang Liu, Yangmuzi Zhang, Dacheng Tao, and Larry S Davis.  Truncated cauchy non-negative matrix factorization.  _IEEE Transactions on pattern analysis and machine intelligence_ , 41(1):246–259, 2017.
  * He and Shao (1996) Xuming He and Qi-Man Shao.  A general bahadur representation of m-estimators and its application to linear regression with nonstochastic designs.  _The Annals of Statistics_ , 24(6):2608–2630, 1996.
  * Huber (1964) Peter J. Huber.  Robust estimation of a location parameter.  _Annals of Mathematical Statistics_ , 35(1):492–518, 1964.
  * Jiao et al. (2021) Yuling Jiao, Guohao Shen, Yuanyuan Lin, and Jian Huang.  Deep nonparametric regression on approximately low-dimensional manifolds.  _arXiv preprint arXiv:2104.06708_ , 2021.
  * Khamaru and Wainwright (2019) Koulik Khamaru and Martin J Wainwright.  Convergence guarantees for a class of non-convex and non-smooth optimization problems.  _Journal of Machine Learning Research_ , 20(154):1–52, 2019.
  * Klusowski et al. (2019) Jason M Klusowski, Dana Yang, and WD Brinda.  Estimating the coefficients of a mixture of two linear regressions by expectation maximization.  _IEEE Transactions on Information Theory_ , 65(6):3515–3524, 2019.
  * Knight (1998) Keith Knight.  Limiting distributions for $l_{1}$ regression estimators under general conditions.  _Annals of Statistics_ , 26(2):755–770, 1998\.
  * Koenker (2005) Roger Koenker.  _Quantile regression_.  Cambridge University Press, New York, 2005.
  * Koenker and Bassett Jr (1978) Roger Koenker and Gilbert Bassett Jr.  Regression quantiles.  _Econometrica_ , pages 33–50, 1978.
  * Lam and Cheng (2021) Clifford Lam and Wenyu Cheng.  Robust mean and eigenvalues regularized covariance matrix estimation.  _London School of Economics and Political Science_ , 2021.
  * Lecué and Mendelson (2013) Guillaume Lecué and Shahar Mendelson.  Learning subgaussian classes: upper and minimax bounds.  _Topics in Learning Theory-Societe Mathematique de France,(S. Boucheron and N. Vayatis Eds.)_ , 2013.
  * Lee et al. (2020) Kyungjae Lee, Hongjun Yang, Sungbin Lim, and Songhwai Oh.  Optimal algorithms for stochastic multi-armed bandits with heavy tailed rewards.  In _Advances in Neural Information Processing Systems_ , volume 33, 2020.
  * Lerasle (2019) Matthieu Lerasle.  Lecture notes: Selected topics on robust statistical learning theory.  _arXiv preprint arXiv:1908.10761_ , 2019.
  * Liu and Tao (2014) Tongliang Liu and Dacheng Tao.  On the robustness and generalization of cauchy regression.  In _2014 4th IEEE International Conference on Information Science and Technology_ , pages 100–105. IEEE, 2014.
  * Loh (2017) Po-Ling Loh.  Statistical consistency and asymptotic normality for high-dimensional robust $m$-estimators.  _The Annals of Statistics_ , 45(2):866–896, 2017\.
  * Lugosi and Mendelson (2019) Gábor Lugosi and Shahar Mendelson.  Mean estimation and regression under heavy-tailed distributions: A survey.  _Foundations of Computational Mathematics_ , 19(5):1145–1190, 2019.
  * Mathieu and Minsker (2021) Timothée Mathieu and Stanislav Minsker.  Excess risk bounds in robust empirical risk minimization.  _Information and Inference: A Journal of the IMA_ , 2021.
  * McCullagh and Nelder (1989) P McCullagh and John A Nelder.  _Generalized Linear Models_ , volume 37.  CRC Press, 1989.
  * Mei et al. (2018) Song Mei, Yu Bai, Andrea Montanari, et al.  The landscape of empirical risk for nonconvex losses.  _Annals of Statistics_ , 46(6A):2747–2774, 2018\.
  * Minsker (2018) Stanislav Minsker.  Sub-gaussian estimators of the mean of a random matrix with heavy-tailed entries.  _The Annals of Statistics_ , 46(6A):2871–2903, 2018.
  * Nemirovskij and Yudin (1983) Arkadij Semenovič Nemirovskij and David Borisovich Yudin.  Problem complexity and method efficiency in optimization.  1983\.
  * Ostrovskii and Bach (2021) Dmitrii M Ostrovskii and Francis Bach.  Finite-sample analysis of $m$-estimators using self-concordance.  _Electronic Journal of Statistics_ , 15(1):326–391, 2021.
  * Resnick (2007) Sidney I Resnick.  _Heavy-tail phenomena: probabilistic and statistical modeling_.  Springer Science & Business Media, 2007.
  * Riphahn et al. (2003) Regina T Riphahn, Achim Wambach, and Andreas Million.  Incentive effects in the demand for health care: a bivariate panel count data estimation.  _Journal of applied econometrics_ , 18(4):387–405, 2003.
  * Schreuder et al. (2020) Nicolas Schreuder, Victor-Emmanuel Brunel, and Arnak Dalalyan.  A nonasymptotic law of iterated logarithm for general m-estimators.  In _International Conference on Artificial Intelligence and Statistics_ , pages 1331–1341. PMLR, 2020.
  * Shao and Zhang (2021) Qi-Man Shao and Zhuo-Song Zhang.  Berry–esseen bounds for multivariate nonlinear statistics with applications to m-estimators and stochastic gradient descent algorithms.  _arXiv preprint arXiv:2102.04923_ , 2021.
  * Shen et al. (2021a) Guohao Shen, Yuling Jiao, Yuanyuan Lin, and Jian Huang.  Non-asymptotic excess risk bounds for classification with deep convolutional neural networks.  _arXiv preprint arXiv:2105.00292_ , 2021a.
  * Shen et al. (2021b) Guohao Shen, Yuling Jiao, Yuanyuan Lin, and Jian Huang.  Robust nonparametric regression with deep neural networks.  _arXiv preprint arXiv:2107.10343_ , 2021b.
  * Sun (2021) Qiang Sun.  Do we need to estimate the variance in robust mean estimation?  _arXiv preprint arXiv:2107.00118_ , 2021.
  * Sun et al. (2020) Qiang Sun, Wen-Xin Zhou, and Jianqing Fan.  Adaptive huber regression.  _Journal of the American Statistical Association_ , 115(529):254–265, 2020.
  * Taheri et al. (2021) Mahsa Taheri, Fang Xie, and Johannes Lederer.  Statistical guarantees for regularized neural networks.  _Neural Networks_ , 142:148–161, 2021.
  * Tukey (1960) John W Tukey.  A survey of sampling from contaminated distributions.  _Contributions to probability and statistics_ , pages 448–485, 1960\.
  * Van der Vaart (1998) Aad W Van der Vaart.  _Asymptotic statistics_ , volume 3.  Cambridge university press, 1998.
  * Vershynin (2018) Roman Vershynin.  _High-dimensional probability: An introduction with applications in data science_ , volume 47.  Cambridge university press, 2018.
  * Wang et al. (2022) Lili Wang, Chao Zheng, Wen Zhou, and Wen-Xin Zhou.  A new principle for tuning-free huber regression.  _Statistica Sinica_ , 32:1–25, 2022.
  * Xu et al. (2020) Yi Xu, Shenghuo Zhu, Sen Yang, Chi Zhang, Rong Jin, and Tianbao Yang.  Learning with non-convex truncated losses by sgd.  In _Uncertainty in Artificial Intelligence_ , pages 701–711. PMLR, 2020.
  * Yang et al. (2021) Xiaowei Yang, Shuang Song, and Huiming Zhang.  Law of iterated logarithm and model selection consistency for generalized linear models with independent and dependent responses.  _Frontiers of Mathematics in China_ , pages 1–32, 2021.
  * Yi et al. (2020) Mingyang Yi, Ruoyu Wang, and Zhi-Ming Ma.  Non-asymptotic analysis of excess risk via empirical risk landscape.  _arXiv preprint arXiv:2012.02456_ , 2020.
  * Zhang and Jia (2022) Huiming Zhang and Jinzhu Jia.  Elastic-net regularized high-dimensional negative binomial regression: Consistency and weak signals detection.  _Statistica Sinica_ , 32:181–207, 2022.
  * Zhang and Zhou (2018) Lijun Zhang and Zhi-Hua Zhou.  $\ell_{1}$-regression with heavy-tailed distributions.  In _NeurIPS_ , 2018.
  * Zhu and Zhou (2021) Ziwei Zhu and Wenjing Zhou.  Taming heavy-tailed features by shrinkage.  In _International Conference on Artificial Intelligence and Statistics_ , pages 3268–3276. PMLR, 2021.

## Appendix

The appendix includes the proofs of the lemmas, corollaries and theorems in
the main body.

##  A Proofs

###  A.1 Useful lemmas and results

The following result is in Corollary 4.2.13 in Vershynin (2018).

######  Lemma 7 (Covering numbers of the Euclidean ball)

For any $\kappa>0$, the covering numbers of the $n$-dimensional unit Euclidean
ball $B_{2}^{n}$ satisfy
$\left(\frac{1}{\kappa}\right)^{n}\leq{N}\left(B_{2}^{n},\kappa\right)\leq\left(\frac{2}{\kappa}+1\right)^{n}.$

###  A.2 The proof of Theorem 1

Proof  By the definition of $\hat{\theta}_{n}$, we have

$\displaystyle{{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\hat{\theta}_{n})=\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi_{\lambda}[\alpha
l({Y_{i}},{{X}}_{i},\hat{\theta}_{n})]\leq\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi_{\lambda}[\alpha
l({Y_{i}},{{X}}_{i},\theta^{*})]=:{{\hat{R}}_{\psi_{\alpha}\circ
l}}({\theta^{*}})$ (28)

which yields

$\displaystyle{R_{l}}(\hat{\theta}_{n})-R_{l}({\theta^{*}})$
$\displaystyle={\left({{R_{l}}(\hat{\theta}_{n})-{{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\hat{\theta}_{n})}\right)}+{\left({{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\hat{\theta}_{n})-{{\hat{R}}_{\psi_{\lambda},l,\alpha}}({\theta^{*}})\right)}+{\left({{{\hat{R}}_{\psi_{\lambda},l,\alpha}}({\theta^{*}})-{R_{l}}({\theta^{*}})}\right)}$
$\displaystyle\leq\left({{{\hat{R}}_{\psi_{\lambda},l,\alpha}}({\theta^{*}})-{R_{l}}({\theta^{*}})}\right)+\left({{R_{l}}(\hat{\theta}_{n})-{{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\hat{\theta}_{n})}\right).$
(29)

Under the two high provability events $\mathcal{E}_{1}(l,\lambda)$ in Lemma 8
and $\mathcal{E}_{2}(l,\lambda)$ in Lemma 9 below, we have by inequality (A.2)

$\displaystyle~{}~{}~{}~{}{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})\leq\left(\frac{f(\alpha)}{\alpha}{R_{\lambda\circ
l}}({\theta^{*}})+\frac{1}{{n\alpha}}\log\frac{1}{\delta}\right)$
$\displaystyle+2\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta}$
$\displaystyle=2\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{{(c_{2}+1)f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta^{2}}.$ (30)

with probability at least $1-2\delta$.

For the last two terms in (A.2), we put
$\frac{{(c_{2}+1)f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)=\frac{1}{{n\alpha}}\log\frac{{N({\Theta},\kappa)}}{{{\delta^{2}}}}$,
i.e. the variance term equals to the bias term. Then we get
$\alpha={f^{-1}}\left({\frac{1}{{n({c_{2}}+1)}}[\mathop{\sup}\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta){]^{-1}}\log\frac{{N(\Theta,\kappa)}}{{{\delta^{2}}}}}\right)$. From
Lemma 7 and (C.0), we have

$\displaystyle\log(\frac{1}{{{\delta^{2}}}})+p\log\left({\frac{r}{\kappa}}\right)$
$\displaystyle\leq\log\frac{{N(\Theta,\kappa)}}{{{\delta^{2}}}}\leq\log\frac{{N({B_{r}},\kappa)}}{{{\delta^{2}}}}$
$\displaystyle\leq\log\left({\frac{1}{{{\delta^{2}}}}{{\left({1+\frac{{2r}}{\kappa}}\right)}^{p}}}\right)=\log(\frac{1}{{{\delta^{2}}}})+p\log\left({1+\frac{{2r}}{\kappa}}\right).$
(31)

By (C.1), the $f(t)$ is increasing function of $t$. We get

${\underline{\alpha}_{n}}:={f^{-1}}\left({\frac{{\log({\delta^{-2}})+p\log\left({r/\kappa}\right)}}{{n({c_{2}}+1){\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)}}}\right)\leq\alpha\leq{f^{-1}}\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2r/\kappa}\right)}}{{n({c_{2}}+1){\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)}}}\right):={{\bar{\alpha}}_{n}}.$

From this upper bound and lower bound of $\alpha$, (A.2) implies

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq
2\kappa{{\rm{E}}}[H(Y,X)]+\frac{{{c_{2}}f({{\bar{\alpha}}_{n}}\kappa)}}{\underline{\alpha}_{n}}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{2}{{n\alpha}}\log\frac{N(\Theta,\kappa)}{\delta^{2}}$
$\displaystyle[\text{By}~{}\eqref{eq:covering}]~{}$ $\displaystyle\leq
2\kappa{{\rm{E}}}[H(Y,X)]+\frac{{{c_{2}}f({{\bar{\alpha}}_{n}}\kappa)}}{\underline{\alpha}_{n}}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{2}{{n{\underline{\alpha}_{n}}}}[{\log({\delta^{-2}})+p\log\left({1+2r/\kappa}\right)}],$

with probability at least $1-2\delta$. Let $\kappa=\frac{1}{n}$, we obtain

$R_{l}(\hat{\theta}_{n})-R_{l}(\theta^{*})\leq\frac{{2{{\rm{E}}}[H(Y,X)]}}{{n}}+\frac{{{c_{2}}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}}}{{{\underline{\alpha}_{n}}}}f(\frac{{\bar{\alpha}}_{n}}{n})+\frac{2}{{n{\underline{\alpha}_{n}}}}[{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}].$

######  Lemma 8 (Generalization/concentration error bound)

For general loss function $l(\cdot,\cdot,\cdot)$, under (C.4), we have

$P\\{\mathcal{E}_{1}(l,\lambda)\\}\geq 1-\delta,$

where
$\mathcal{E}_{1}(l,\lambda):=\\{{{{\hat{R}}_{\psi_{\lambda},l,\alpha}}({\theta^{*}})-{{R}_{l}}({\theta^{*}})\leq\frac{f(\alpha)}{\alpha}{R_{\lambda\circ
l}}({\theta^{*}})+\frac{1}{{n\alpha}}\log\frac{1}{\delta}}\\}$.

Proof  The proof is based on bounding the exponential moment of
$n\alpha{{\hat{R}}_{\psi_{\lambda},l,\alpha}}({\theta^{*}})$ from Markov’s
inequality. Applying the upper truncated function in (4), $1+x\leq e^{x}$, one
has

$\displaystyle{\rm{E}}{e^{n\alpha{{\hat{R}}_{\psi_{\alpha}\circ
l}}({\theta^{*}})}}$
$\displaystyle={\rm{E}}{e^{\sum\limits_{i=1}^{n}\psi[\alpha
l({Y_{i}},X_{i},{\theta^{*}})]}}\leq{\rm{E}}\\{\prod\limits_{i=1}^{n}{[1+\alpha
l({Y_{i}},X_{i},{\theta^{*}})+\lambda[\alpha
l({Y_{i}},X_{i},\theta^{*})]]}\\}$ $\displaystyle[\text{By independence}]~{}$
$\displaystyle=\prod\limits_{i=1}^{n}{\\{{\rm{E}}[1+\alpha
l({Y}_{i},X_{i},{\theta^{*}})+\lambda[\alpha
l({Y_{i}},X_{i},\theta^{*})]]\\}}$ $\displaystyle\leq
e^{\alpha\sum_{i=1}^{n}{\rm{E}}l(Y_{i},X_{i},{\theta^{*}})+\sum_{i=1}^{n}{\rm{E}}\\{\lambda[\alpha
l({Y_{i}},X_{i},\theta^{*})]\\}}\leq{e^{n[\alpha{R_{l}}({\theta^{*}})+f(\alpha){R_{\lambda\circ
l}}({\theta^{*}})]}}.$

Via Markov’s inequality with the exponential transform, it gives

$\displaystyle~{}~{}~{}~{}P\\{\mathcal{E}_{1}^{c}(l,c)\\}=P\\{{{\hat{R}}_{\psi_{\lambda},l,\alpha}}({\theta^{*}})>[{R_{l}}({\theta^{*}})+\frac{f(\alpha)}{\alpha}{R_{\lambda\circ
l}}({\theta^{*}})]+\frac{{\log(1/\delta)}}{{n\alpha}}\\}$
$\displaystyle=P\\{{e^{n\alpha{{\hat{R}}_{\psi_{\lambda},l,\alpha}}({\theta^{*}})}}>{e^{n\alpha[{R_{l}}({\theta^{*}})+\frac{f(\alpha)}{\alpha}{R_{\lambda\circ
l}}({\theta^{*}})]+\log(1/\delta)}}\\}\leq\frac{{{\rm{E}}{e^{n\alpha{{\hat{R}}_{\psi_{\lambda},l,\alpha}}({\theta^{*}})}}}}{{{e^{n\alpha[{R_{l}}({\theta^{*}})+\frac{f(\alpha)}{\alpha}{R_{\lambda\circ
l}}({\theta^{*}})]+\log(1/\delta)}}}}\leq\delta.$

######  Lemma 9 (Generalization error bound)

For any $\kappa>0$, under (C.1)-(C.4), one has

$P\\{\mathcal{E}_{2}(l,\lambda)\\}\geq 1-\delta,$

where
$\mathcal{E}_{2}(l,\lambda):=\\{{{R}_{l}}(\hat{\theta}_{n})-{{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\hat{\theta}_{n})\leq
2\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta}\\}$.

Proof  Let $\mathcal{N}({\Theta},\kappa)$ be an $\kappa$-net of ${\Theta}$ and
denote its covering number as ${N}({\Theta},\varepsilon)$. For each
$\hat{\theta}_{n}\in{\Theta}$, the definition of $\kappa$-net implies that
there exists a $\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)$ satisfying
$\|\hat{\theta}_{n}-\tilde{\theta}\|_{{2}}\leq\kappa$, by locally Lipschitz
condition (C.2), we obtain

$\displaystyle l({Y_{i}},{{X}}_{i},\hat{\theta}_{n})\geq
l({Y_{i}},{{X}}_{i},\tilde{\theta})-\kappa
H({Y_{i}},{{X}}_{i}),~{}i=1,2,\cdots,n.$

Since $\psi_{\lambda}(\cdot)$ is non-decreasing and the last inequality gives

$\displaystyle{{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\hat{\theta}_{n})$
$\displaystyle=\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi_{\lambda}[\alpha
l({Y_{i}},X_{i},\hat{\theta}_{n})]\geq\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi_{\lambda}[\alpha
l({Y_{i}},X_{i},\tilde{\theta})-\kappa\alpha H({Y_{i}},{{X}}_{i})].$ (32)

In below, we continue to derive the lower bound in (32) by applying covering
number techniques. From the lower bound in (4), we could estimate the
exponential moment bound:

$\displaystyle~{}~{}~{}~{}{\rm{E}}{e^{-\sum_{i=1}^{n}\psi_{\lambda}[\alpha
l({Y_{i}},X_{i},\tilde{\theta})-\alpha\kappa H({Y_{i}},{{X}}_{i})]}}$
$\displaystyle\leq{\rm{E}}\prod\limits_{i=1}^{n}\\{1-\alpha{\rm{E}}[l({Y_{i}},X_{i},\tilde{\theta})]+\alpha\kappa
H({Y_{i}},{{X}}_{i})+\lambda[\alpha(l({Y_{i}},X_{i},\tilde{\theta})-\kappa
H({Y_{i}},{{X}}_{i}))]\\}$
$\displaystyle\leq\prod\limits_{i=1}^{n}\left\\{1-\alpha{\rm{E}}[l({Y_{i}},X_{i},\tilde{\theta})]+\alpha\kappa{\rm{E}}[H({Y_{i}},{{X}}_{i})]+{\rm{E}}\\{\lambda[\alpha(l({Y_{i}},X_{i},\tilde{\theta})-\kappa
H({Y_{i}},{{X}}_{i}))]\\}\right\\}$ $\displaystyle[\text{(C.1.2)}]$
$\displaystyle\leq\prod\limits_{i=1}^{n}\left\\{1-\alpha{\rm{E}}[l({Y_{i}},X_{i},\tilde{\theta})]+\alpha\kappa{\rm{E}}[H({Y_{i}},{{X}}_{i})]+{c_{2}}{\rm{E}}\\{\lambda[\alpha(l({Y_{i}},X_{i},\tilde{\theta}))]+\lambda[\alpha\kappa
H({Y_{i}},{X_{i}})]\\}\right\\}$ $\displaystyle[\text{(C.1.1)}]$
$\displaystyle\leq{e^{n\alpha\\{{-R_{l}(\tilde{\theta})+\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha)}}{\alpha}R_{\lambda\circ
l}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}}\\}}},$

where the last inequality stems from $1+x\leq e^{x}$. By the last exponential
moment bound, Markov’s inequality shows that, for a fixed
$\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)$,

$\displaystyle P$
$\displaystyle\left\\{\frac{{-1}}{{n\alpha}}\sum\limits_{i=1}^{n}\psi_{\lambda}(\alpha
l({Y_{i}},X_{i},\tilde{\theta})-\kappa\alpha
H({Y_{i}},{{X}}_{i}))>-{R_{l}}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha)}}{\alpha}R_{\lambda\circ
l}(\tilde{\theta})\right.$
$\displaystyle\left.+\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+{\frac{{\log(1/s)}}{{n\alpha}}}\right\\}$
$\displaystyle\leq\frac{{e^{n\alpha\\{{-R_{l}(\tilde{\theta})+\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha)}}{\alpha}R_{\lambda\circ
l}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}}\\}}}}{{{e^{n\alpha[-R_{l}(\tilde{\theta})+\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha)}}{\alpha}R_{\lambda\circ
l}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}]+\log(1/s)}}}}=s\in(0,1).$
(33)

The set $\mathcal{N}({\Theta},\kappa)$ has $N({\Theta},\kappa)$ elements. From
the single bound (A.2) and putting $s=\delta/N({\Theta},\kappa)$, we have for
all $\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)$

$\displaystyle P$
$\displaystyle\left(\bigcup\limits_{\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)}\left\\{\frac{{-1}}{{n\alpha}}\sum\limits_{i=1}^{n}\psi_{\lambda}(\alpha
l({Y_{i}},X_{i},\tilde{\theta})-\kappa\alpha
H({Y},{{X}}))>-{R_{l}}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha)}}{\alpha}R_{\lambda\circ
l}(\tilde{\theta})\right.\right.$
$\displaystyle\left.\left.+\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+{\frac{{\log(1/s)}}{{n\alpha}}}\right\\}\right)$
$\displaystyle\leq N({\Theta},\kappa)\cdot
P\left\\{\frac{{-1}}{{n\alpha}}\sum\limits_{i=1}^{n}\psi_{\lambda}(\alpha
l({Y_{i}},X_{i},\tilde{\theta})-\kappa\alpha
H({Y},{{X}}))\geq-{R_{l}}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha)}}{\alpha}R_{\lambda\circ
l}(\tilde{\theta})\right.$
$\displaystyle\left.+\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+{\frac{{\log(1/s)}}{{n\alpha}}}\right\\}\leq
N({\Theta},\kappa)s=:\delta.$ (34)

Then the complementary set in (A.2) hold with probability at least $1-\delta$.
Inequalities (32) and (A.2) give the following lower bound for all
$\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)$

$\displaystyle~{}~{}~{}~{}{{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\hat{\theta})\geq\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi_{\lambda}[\alpha
l({Y_{i}},X_{i},\tilde{\theta})-\kappa\alpha H({Y},{{X}})]$
$\displaystyle\geq{{R}_{l}}(\tilde{\theta})-\left\\{\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha)}}{\alpha}R_{\lambda\circ
l}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta}\right\\}$
$\displaystyle\geq{{R}_{l}}(\tilde{\theta})-\left\\{\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta}\right\\}.$
(35)

with probability at least $1-\delta$.

It remains to find the lower bound for ${{R}_{l}}(\tilde{\theta})$ in (A.2) by
the error bound of $|{R_{l}}(\hat{\theta}_{n})-{R_{l}}(\tilde{\theta})|$. To
this end, the (C.2) implies

${R_{l}}(\hat{\theta}_{n})-{R_{l}}(\tilde{\theta})\leq{\rm{E}}[H(X,Y)\|\hat{\theta}_{n}-\tilde{\theta}\|_{2}]\leq\kappa{\rm{E}}[H(X,Y)]~{}~{}\hat{\theta},\tilde{\theta}\in{B_{r}},$

which gives
${R_{l}^{n}}(\tilde{\theta})\geq{R_{l}}(\hat{\theta}_{n})-\kappa{\rm{E}}[H(X,Y)]$.
Thus, (A.2) has a further lower bound:

$\displaystyle~{}~{}~{}~{}{{\hat{R}}_{\psi_{\lambda},l,\alpha}}(\hat{\theta}_{n})$
$\displaystyle\geq{{R}_{l}}(\hat{\theta}_{n})-\\{2\kappa{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\tilde{\theta})+\frac{{{c_{2}}f(\alpha\kappa)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta}\\}$

with probability at least $1-\delta$. Then, we conclude Lemma 9.

###  A.3 The proof of Theorem 2

Proof  Let $\kappa=1/n$ in (A.2), we get with probability at least $1-2\delta$

$\displaystyle~{}~{}~{}~{}{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{2}{n}{{\rm{E}}}[H({Y},{{X}})]+\frac{{{c_{2}}f(\alpha/n)}}{\alpha}{{\rm{E}}}\\{\lambda[H(Y,X)]\\}+\frac{{(c_{2}+1)f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)+\frac{1}{{n\alpha}}\log\frac{N(\Theta,1/n)}{\delta^{2}}.$

The $\lambda(x)=|x|^{\beta}/{\beta},~{}{\beta}\in(1,2]$ satisfies weak
triangle inequality and homogeneous inequality

$\displaystyle|x+y|^{\beta}/{\beta}\leq
2^{{\beta}-1}[|x|^{\beta}/{\beta}+|y|^{\beta}/{\beta}],~{}~{}|tx|^{\beta}/{\beta}\leq|t|^{\beta}\cdot|x|^{\beta}/{\beta}$
(36)

by $|a+b|^{\beta}\leq 2^{{\beta}-1}(|a|^{\beta}+|b|^{\beta})$ for ${\beta}>1$.
Thus, we have $c_{2}=2^{{\beta}-1}$ and $f(t)=|t|^{\beta}$ in (C.1).

Using the variance-bias tradeoff, i.e. the variance term equals to the bias
term, put
$\frac{{(c_{2}+1)f(\alpha)}}{\alpha}{\sup}_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)=\frac{1}{{n\alpha}}\log\frac{{N({\Theta},\kappa)}}{{{\delta^{2}}}}$.
Note that $f^{-1}(t)=t^{1/\beta}$ for $t>0$. So we have

$\alpha=f^{-1}\left(\frac{1}{n({c_{2}}+1)}[\sup\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]^{-1}\log\frac{{N(\Theta,\kappa)}}{{\delta^{2}}}\right)=\frac{1}{n^{1/\beta}}\left(\frac{\log[{{N(\Theta,\kappa)}}/{{\delta^{2}}}]}{(2^{{\beta}-1}+1)\sup\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)}\right)^{1/\beta}$.

Observe that $f(t)/t=t^{\beta-1}$ for $t>0$. Then

$\displaystyle~{}~{}~{}~{}\frac{{(c_{2}+1)f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)+\frac{1}{{n\alpha}}\log\frac{N(\Theta,1/n)}{\delta^{2}}=\frac{{2(c_{2}+1)f(\alpha)}}{\alpha}\mathop{\sup}\limits_{\theta\in{\Theta}}R_{\lambda\circ
l}(\theta)$
$\displaystyle=\frac{{2({2^{\beta-1}}+1)}}{{{n^{(\beta-{\rm{1)}}/\beta}}}}{\left({\frac{{\log[N(\Theta,1/n)/{\delta^{2}}]}}{{({2^{\beta-1}}+1)\mathop{\sup}\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)}}}\right)^{(\beta-{\rm{1)}}/\beta}}\mathop{\sup}\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)$
$\displaystyle=\frac{{2({2^{\beta-1}}+1)}}{{{n^{(\beta-{\rm{1)}}/\beta}}}}{\left({\frac{{\log[N(\Theta,1/n)/{\delta^{2}}]}}{{{2^{\beta-1}}+1}}}\right)^{(\beta-{\rm{1)}}/\beta}}{[\mathop{\sup}\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]^{\beta^{-1}}}$
$\displaystyle\leq\frac{{2({2^{\beta-1}}+1)}}{{{n^{(\beta-{\rm{1)}}/\beta}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2r/\kappa}\right)}}{{{2^{\beta-1}}+1}}}\right)^{(\beta-{\rm{1)}}/\beta}}{[\mathop{\sup}\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]^{\beta^{-1}}},$

where the last inequality is by (A.2).

Thus, $\frac{{f(\alpha/n)}}{\alpha}={\alpha}^{\beta-1}\cdot n^{-\beta}$ and
$\lambda(x)=|x|^{\beta}/{\beta}$ show

$\displaystyle~{}~{}~{}~{}\frac{{{c_{2}}f(\alpha/n)}}{\alpha}{{\rm{E}}}\\{\lambda[{H}(Y,X)]\\}=\frac{{{2^{\beta-1}/\beta}}}{{{n^{1-{\beta^{-1}}+\beta}}}}{\left({\frac{{\log[N(\Theta,1/n)/{\delta^{2}}]}}{{({2^{\beta-1}}+1)\mathop{\sup}\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)}}}\right)^{(\beta-{\rm{1)}}/\beta}}{{\rm{E}}}\\{|{H}(Y,X){|^{\beta}}\\}$
$\displaystyle\leq\frac{{{2^{\beta-1}/\beta}}}{{{n^{1-{\beta^{-1}}+\beta}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2r/\kappa}\right)}}{{{2^{\beta-1}}+1}}}\right)^{(\beta-{\rm{1)}}/\beta}}\frac{{{{\rm{E}}}\\{|{H}(Y,X){|^{\beta}}\\}}}{{{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{(\beta-{\rm{1)}}/\beta}}}},$

where the last inequality is by (A.2). Then we have

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{2}{n}{{\rm{E}}}[H({Y},{{X}})]+\frac{1}{{{n^{{{\frac{\beta-1}{\beta}}}}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2r/\kappa}\right)}}{{{2^{\beta-1}}+1}}}\right)^{{\textstyle{{\beta-{\rm{1}}}\over\beta}}}}$
$\displaystyle\cdot\left[{\frac{{{2^{\beta-1}}{{\rm{E}}}|{H}(Y,X)|^{\beta}/\beta}}{{{n^{\beta}}{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{{\textstyle{{\beta-{\rm{1}}}\over\beta}}}}}}+2({2^{\beta-1}}+1){{[\mathop{\sup}\limits_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{{\beta^{-1}}}}}\right],$

with probability at least $1-2\delta$. We get the desired result by letting
$\beta={1+\varepsilon}$.

###  A.4 More results, remarks and proofs for Section 3

####  A.4.1 Quantile regression with heterogeneous regression noises

Similar to Corollary 4 with second moment assumption, assume the sample
$\\{(X_{i},Y_{i})\\}_{i=1}^{n}$ satisfies that $\\{X_{i}\\}_{i=1}^{n}$ are
i.i.d. while $\\{Y_{i}\\}_{i=1}^{n}$ are non-identically distributed obeying
(11) with heterogeneous regression noises (Koenker, 2005). The sample-
dependent parameter $\theta_{n}^{*}$ is defined as the minimizer

$\theta_{n}^{*}:=\mathop{\arg\min}\limits_{\theta\in{\Theta}}\mathrm{E}[\frac{1}{n}\sum\limits_{i=1}^{n}{{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\theta)}]$
and
${{R}_{\rho_{\tau}}^{n}}(\theta):=\mathrm{E}[\frac{1}{n}\sum_{i=1}^{n}{{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\theta)}]$
for ${\theta\in{\Theta}}$.

To guarantee high probability $\ell_{2}$-estimation error, it requires the
minimal eigenvalue condition of the covariance matrix for the design and the
positivity of the density functions of heterogeneous regression noises
$\\{{f_{\varepsilon_{i}}}(t)\\}_{i=1}^{n}$.

###### Corollary 10

Let
$\sigma_{R}^{\rm{2}}(n):={\sup}_{\theta\in{\Theta}}\frac{1}{n}\sum_{i=1}^{n}\mathrm{E}[{\rho^{2}_{\tau}}({Y}_{i}-{{X}_{i}}^{\top}\theta)]<\infty$
and ${{\rm{E}}}\left\|{{X}}\right\|_{{2}}^{2}<\infty$. With probability at
least $1-2\delta$,one has

${R_{{\rho_{\tau}}}^{n}}(\hat{\theta}_{n})-{R_{{\rho_{\tau}}}^{n}}({\theta_{n}^{*}})\leq\frac{{2{l_{\tau}}}}{n}{\rm{E}}{\left\|{{X}}\right\|_{{2}}}+\sqrt{\frac{{2\log({\delta^{-2}})+2p\log\left({1+2rn}\right)}}{n}}\left[\sqrt{3}\sigma_{R}(n)+{\frac{{{l_{\tau}^{2}}{\rm{E}}\left\|{{X}}\right\|_{{2}}^{2}}}{{\sqrt{3}\sigma_{R}(n)n^{2}}}}\right].$

Further, consider: (Qa) The $\\{{\varepsilon}_{i}\\}_{i=1}^{n}$ are the
continuous r.v.s and its conditional density functions
$\mathop{\inf}\limits_{|t|<\infty}\min\limits_{1\leq i\leq
n}{f_{\varepsilon_{i}}}(t)>{c_{f}},$ where ${c_{f}}$ is a positive constant;
(Qb) Assume that there exists a positive constant ${c_{\Sigma}}$ such that
${\rm{E}}[XX^{\top}]\succ\mathrm{I}_{p}{c_{\Sigma}}.$ Then, we have with
probability at least $1-2\delta$

${\|{\hat{\theta}_{n}-\theta_{n}^{*}}\|_{2}^{2}}\leq\frac{1}{{{c_{f}}{c_{\Sigma}}}}\left\\{\frac{{2{l_{\tau}}}}{n}{\rm{E}}{\left\|{{X}}\right\|_{{2}}}+\sqrt{\frac{{2\log({\delta^{-2}})+2p\log\left({1+2rn}\right)}}{n}}\left[\sqrt{3}\sigma_{R}(n)+{\frac{{{l_{\tau}^{2}}{\rm{E}}\left\|{{X}}\right\|_{{2}}^{2}}}{{\sqrt{3}\sigma_{R}(n)n^{2}}}}\right]\right\\}.$

###  A.5 Remarks and proofs for GLMs

###### Remark 11

The (14) is originally derived from negative log-likelihood functions of
exponential family. The exponential family contains many sub-exponential and
sub-Gaussian distributions such as binomial, Poisson, negative binomial,
Normal, Gamma distributions (McCullagh and Nelder, 1989). Let $\nu(\cdot)$ be
some dominated measure and $b(\cdot)$ be a function. Consider a $Y$ follows
the distribution of the natural exponential families $P_{\eta}$ indexed by
parameter ${\eta}$.

$P_{\eta}(dy)=c(y)\exp\\{y{\eta}-b({\eta})\\}\nu(dy),$ (37)

where the function $c(y)>0$ is free of ${\eta}\in\Xi:=\\{{\eta}:\smallint
c(y)\exp\\{y{\eta}\\}\nu(dy)<\infty\\}.$

###### Remark 12

For GLMs, the link function $u(x)$ is canonical, i.e. $u(x)=x$, whence we can
choose ${g_{A}}(\cdot)\equiv 1$ in (G1). Note that in this case $k(t)=b(t)$, a
choice of ${h_{A}}(\cdot)$ in condition (G2) is derived by $\ddot{b}(t)>0$

$\dot{k}(x^{\top}\theta)\leq\dot{b}(r\left\|{x}\right\|_{2})=:{h_{r}}(x),~{}\text{for}~{}\|\theta\|_{2}\leq
r~{}\text{due to}~{}|{x^{\top}\theta}|\leq r\left\|{x}\right\|_{2}.$

For GLM with non-canonical link function $u(x)$, we first choose
${g_{A}}(\cdot)$ in condition (G1) by
$\dot{u}(x^{\top}\theta)\leq\dot{u}(r\left\|x\right\|_{2})=:{g_{r}}(x),~{}\text{for}~{}\|\theta\|_{2}\leq
r$ due to $|{x^{\top}\theta}|\leq r\left\|x\right\|_{2}.$ Under the (C.0) and
(G.1), it implies (G.2) with $A=r$ and
${h_{r}}(x)={g_{r}}(x)\dot{b}(u(r\left\|x\right\|_{2}))$ by the following
inequality

$\dot{k}(x^{\top}\theta)=\dot{u}(x^{\top}\theta)\dot{b}(u(x^{\top}\theta))\leq{g_{r}}(x)\dot{b}(u(r\left\|x\right\|_{2})):={h_{r}}(x),~{}\text{for}~{}\|\theta\|_{2}\leq
r.$

Suppose the input $\\{X_{i}\\}_{i=1}^{n}$ is i.i.d. drawn from $X$ and $X$ is
bounded (see Yang et al. (2021)). Under the (C.0), the $k(X^{\top}\theta)$ and
$u(X^{\top}\theta)$ are also bounded, then (G.4) is true under the finite
second moments of output

$R_{{l^{2}}}(\theta)={{\rm{E}}}{[k(X^{\top}\theta)-Yu(X^{\top}\theta)]^{2}}\leq
2{{\rm{E}}}{[k(X^{\top}\theta)]^{2}}+2{{\rm{E}}}{[Yu(X^{\top}\theta)]^{2}}\leq{C_{1}}+{C_{2}}{{\rm{E}}}Y^{2}<\infty$

for $\|\theta\|_{2}\leq r$, where $C_{1}$ and $C_{2}$ are some positive
constants.

From (37), one can formally derive the quasi-GLMs loss in (14). Indeed, given
$\\{X_{i}\\}_{i=1}^{n}$, the conditional likelihood function of
$\\{Y_{i}\\}_{i=1}^{n}$ is the product of $n$ terms in (37) with
$\eta_{i}:=u(X_{i}^{\top}\theta)$, and the average negative log-likelihood
function is

${{\hat{R}}_{l}^{n}}(\theta):=\frac{-1}{n}\sum_{i=1}^{n}[Y_{i}u({X}_{i}^{\top}\theta)-b(u({X}_{i}^{\top}\theta))]=\frac{1}{n}\sum_{i=1}^{n}l({Y_{i}},X_{i}^{\top}\theta),~{}\theta\in\mathbb{R}^{p}$.

Fix a $\eta\in{\Theta}$, we compute the gradient for the loss function in (14)

$\displaystyle{\dot{l}(y,x,\eta)}:=\nabla_{\eta}l(y,x^{\top}\eta)=[-y\dot{u}(x^{\top}\eta)+\dot{k}(x^{\top}\eta)]x^{\top}.$
(38)

From (9), ${H}(y,x)$ in Theorem 2 is given by

$\displaystyle{\sup}_{{\eta_{1}},{\eta_{2}}\in\Theta}\|{\dot{l}[y,x,(t{\eta_{2}}+(1-t){\eta_{1}})]}\|_{2}$
$\displaystyle\leq\mathop{\sup}\limits_{\|\theta\|_{2}\leq
A}|-y\dot{u}(x^{\top}\eta)+\dot{k}(x^{\top}\eta)|\cdot\left\|x\right\|_{2}$
$\displaystyle\leq[|y|{g_{A}}(x)+{h_{A}}(x)]\left\|x\right\|_{2}:={H}(y,x)$
(39)

which implies the excess risk bound in Corollary 5.

Robust logistic regression. We have $u(t)=t,~{}k(t)=\log(1+e^{t})$ and
$y\in\\{0,1\\}.$ Note that ${H}(y,x)$ in Theorem 2 is given by

${\sup}_{{\eta_{1}},{\eta_{2}}\in\Theta}\|{\dot{l}[y,x,(t{\eta_{2}}+(1-t){\eta_{1}})]}\|_{2}\leq\mathop{\sup}\limits_{s}|y+\dot{k}(s)|\cdot\left\|x\right\|_{2}\leq
2\left\|x\right\|_{2}:={H}(y,x)$

from (9) and (A.5), we have
$|{H}(y,x){|^{1+\varepsilon}}={2^{1+\varepsilon}}{\left\|x\right\|_{2}^{1+\varepsilon}}$.
Thus Theorem 2 gives

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{{4}{{\rm{E}}}{\left\|{X}\right\|_{2}}}{n}+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}\cdot$
$\displaystyle\left[\frac{2({2^{\varepsilon}}+1)}{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{-{(\varepsilon+1)^{-1}}}}+\frac{{2^{1+2\varepsilon}}{{\rm{E}}}{\left\|{X}\right\|_{2}^{1+\varepsilon}/{({1+\varepsilon})}}}{{n^{1+\varepsilon}{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right].$

Robust negative binomial regression. The connection of $u(\cdot)$ and
$k(\cdot)$ of NBR is $u(t)=t-\log(\eta+e^{t})$ and
$k(t)=\eta\log(\eta+e^{t})$. From (9) and (A.5), we have via Theorem 2

$\displaystyle{\sup}_{{\eta_{1}},{\eta_{2}}\in\Theta}\|{\dot{l}[y,x,(t{\eta_{2}}+(1-t){\eta_{1}})]}\|_{2}$
$\displaystyle\leq\mathop{\sup}\limits_{s}|y+\dot{k}(s)|\cdot\left\|x\right\|_{2}$
$\displaystyle=(y+\eta)\left\|x\right\|_{2}:={H}(y,x),~{}y\geq 0.$ (40)

Hence, we have

$|{H}(y,x){|^{1+\varepsilon}}=|(y+\eta)\left\|x\right\|_{2}|^{1+\varepsilon}\leq{2^{\varepsilon}}[{\left\|yx\right\|_{2}^{1+\varepsilon}}+{(\eta\left\|x\right\|_{2})^{1+\varepsilon}}],$

and thus Theorem 2 gives

$\displaystyle{R_{l}}(\hat{\theta}_{n})-{R_{l}}({\theta^{*}})$
$\displaystyle\leq\frac{2{{\rm{E}}}[(Y+\eta)\left\|X\right\|_{2}]}{n}+\frac{1}{{{n^{\varepsilon/(\varepsilon+1)}}}}{\left({\frac{{\log({\delta^{-2}})+p\log\left({1+2rn}\right)}}{{{2^{\varepsilon}}+1}}}\right)^{\varepsilon/(\varepsilon+1)}}\cdot$
$\displaystyle\left[\frac{2({2^{\varepsilon}}+1)}{{[{\sup}_{\theta\in\Theta}R_{\lambda\circ
l}(\theta)]}^{-{(\varepsilon+1)^{-1}}}}+\frac{{{{{4^{\varepsilon}}}{{\rm{E}}}[{\left\|YX\right\|_{2}^{1+\varepsilon}}+{(\eta\left\|X\right\|_{2})^{1+\varepsilon}/{({1+\varepsilon})}}]}}}{n^{1+\varepsilon}{{{[{\sup_{\theta\in\Theta}}R_{\lambda\circ
l}(\theta)]}^{\varepsilon/(\varepsilon+1)}}}}\right].$

###  A.6 The proof of Theorem 6 and related examples

In (16), let
$\|\theta\|_{\mathrm{F}}:=\sqrt{\sum_{l=0}^{L}\left\|W^{l}\right\|_{\mathrm{F}}^{2}}$
with
$\left\|W^{l}\right\|_{\mathrm{F}}:=\sqrt{\sum_{k=1}^{N_{l+1}}\sum_{j=1}^{N_{l}}(W_{kj}^{l})^{2}}$.
One can check Lipschitz property of DNNs
$\left|{{f_{{\theta_{1}}}}(x)-{f_{{\theta_{2}}}}(x)}\right|\leq
c_{\mathrm{Lip}}({x})\left\|{{\theta_{2}}-{\theta_{1}}}\right\|_{\mathrm{F}}$
for a given function $c_{\mathrm{Lip}}({x})$ by following lemma.

######  Lemma 13 (Lipschitz property of DNNs, Proposition 6 in Taheri et al.
(2021))

Assume that the activation functions
$\sigma_{l}:\mathbb{R}^{N_{l}}\rightarrow\mathbb{R}^{N_{l}}$ are $a$-Lipschitz
with respect to the Euclidean norms on their input and output spaces for
$l=0,1,\dots,L$. Then for every ${x}\in\mathbb{R}^{p}$ and parameters
$\theta=\left(W_{0},\ldots,W_{L}\right),\gamma=\left(V_{0},\ldots,V_{L}\right)$
in (16), it holds that

$\left|f_{\theta}(x)-f_{\gamma}(x)\right|\leq
c_{\mathrm{Lip}}({x})\|\theta-\gamma\|_{\mathrm{F}}$

with
$c_{\mathrm{Lip}}({x}):=2a^{L}\sqrt{L}\|{x}\|_{2}\max_{l\in\\{0,\ldots,L\\}}\prod_{j\in\\{0,\ldots,L\\},j\neq
l}\sigma_{\max}(W^{j})\vee\sigma_{\max}(V^{j}).$

Motivated by Lemma 13 with Lipschitz function $c_{\mathrm{Lip}}({x})$, we
assume that the parameter set has $F$-norm ball with radius $K$ and maximum
spectral norm bounded by $M$,

$\Theta:=\\{\theta:=(W_{1},\ldots,W_{L}):\|\theta\|_{\mathrm{F}}\leq
K~{}\text{and}~{}\max_{0\leq l\leq
L}{\sigma_{\max}}({W_{l}})\leq{M}\\}\subseteq\mathbb{R}^{\sum_{l=0}^{L}N_{l+1}N_{l}}.$
(41)

Suppose that $l(\cdot,\cdot)$ satisfies Lipschitz condition with a function
$D(x,y)$:

$|l(y,f_{\theta_{2}}(x))-l(y,f_{\theta_{1}}(x))|\leq
D(x,y){{|f_{\theta_{2}}(x)-f_{\theta_{1}}(x)|}},~{}\theta_{1},~{}\theta_{2}\in\Theta$.

From Lemma 13, the loss function has Lipschitz property

$|l(y,f_{\theta_{2}}(x))-l(y,f_{\theta_{1}}(x))|\leq
D(x,y)c_{\mathrm{Lip}}({x})\|{\theta_{2}}-{\theta_{1}}\|_{\mathrm{F}}\leq
2(aM)^{L}\sqrt{L}D(x,y)\|{x}\|_{2}\|{\theta_{2}}-{\theta_{1}}\|_{\mathrm{F}}$,

which shows that

${H(y,x):=2(aM)^{L}\sqrt{L}\|{x}\|_{2}D(x,y)~{}\text{in}~{}(C.1).}$ (42)

We conclude Theorem 6 by using Theorem 2.

For robust DNN LAD regression, we have $D(x,y)=1$ and thus
$H(y,x)=2(aM)^{L}\sqrt{L}\|{x}\|_{2}$.

For robust DNN logistic regression, it gives $D(x,y)=y+1\leq 2$ and
$H(y,x)=4(aM)^{L}\sqrt{L}\|{x}\|_{2}$.

For robust DNN NBR, we get $D(x,y)=y+\eta$ and
$H(y,x)=2(aM)^{L}\sqrt{L}\|{x}\|_{2}(y+\eta)$.

###  A.7 The proof of Corollary 10

The proof of Corollary 10 is based on following Lemmas 14. The proof of Lemma
14 is grounded on the two high provability events ${E_{1}}$ in Lemma 15 and
${E_{2}}$ in Lemma 16. The proof of Lemma 15 is also the same as Lemma 8, and
we omit their proofs since the arguments are the same as Theorem 1.

######  Lemma 14 (Basic excess risk bound)

Let $\tau\in(0,1)$ and $\kappa>0$. Define $l_{\tau}:=\max\\{1+\tau,2-\tau\\}$.
With probability at least $1-2\delta$, under (C.0), (Q.1) and (Q.2), we have

${R_{{\rho_{\tau}}}^{n}}(\hat{\theta}_{n})-{R_{{\rho_{\tau}}}^{n}}({\theta_{n}^{*}})\leq
2\kappa l_{\tau}{\rm{E}}{\left\|{{X}}\right\|_{{2}}}+\alpha{(\kappa
l_{\tau})^{2}}{\rm{E}}\left\|{{X}}\right\|_{{2}}^{2}+\frac{3\alpha}{2}\mathop{\sup}\limits_{\theta\in{\Theta}}{R_{\rho_{\tau}^{2}}^{n}}(\theta)+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta^{2}}.$
(43)

######  Lemma 15 (Generalization/concentration error bound)

Under the (Q.3), we have

$P({E_{1}}):=P\left({{{\hat{R}}_{\psi,{\rho_{\tau}},\alpha}}({\theta_{n}^{*}})-{{R}_{{\rho_{\tau}}}^{n}}({\theta_{n}^{*}})\leq\frac{\alpha}{2}{{R}_{\rho_{\tau}^{2}}^{n}}({\theta_{n}^{*}})+\frac{1}{{n\alpha}}\log\frac{1}{\delta}}\right)\geq
1-\delta.$

Next, the following lemma is crucial to complete the proof of Corollary 10.

######  Lemma 16 (Generalization error bound)

For any $\kappa>0$, with (C.0), (Q.1) and (Q.2), the following event (denote
as $E_{2}$) probability exceeding $1-\delta$

${{R}_{{\rho_{\tau}}}^{n}}(\hat{\theta}_{n})-{{\hat{R}}_{\psi,{\rho_{\tau}},\alpha}}(\hat{\theta}_{n})\leq
2\kappa
l_{\tau}{\rm{E}}{\left\|{{X}}\right\|_{{2}}}+\alpha\mathop{\sup}\limits_{\theta\in{\Theta}}{R_{\rho_{\tau}^{2}}^{n}}(\theta)+\alpha{(\kappa
l_{\tau})^{2}}{\rm{E}}\left\|{{X}}\right\|_{{2}}^{2}+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta}.$

Proof  Since $\mathcal{N}({\Theta},\kappa)$ is an $\kappa$-net of ${\Theta}$
with minimal cardinality ${N}({\Theta},\kappa)$ and
$\hat{\theta}_{n}\in{\Theta}$. The definition of $\kappa$-net implies that
there must be a $\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)$ satisfying
$\|\hat{\theta}_{n}-\tilde{\theta}\|_{{2}}\leq\kappa.$ Hence, Cauchy’s
inequality gives

$|{{X}}_{i}^{\top}(\tilde{\theta}-\hat{\theta}_{n})|\leq\kappa\|{{X}}_{i}\|_{{2}}.$
(44)

To prove Lemma 16, we use Knight’s identity (Knight, 1998),

$\rho_{\tau}(u-v)-\rho_{\tau}(u)=-v[\tau-{\rm{1}}(u<0)]+\int_{0}^{v}[{\rm{1}}(u\leq
s)-{\rm{1}}(u\leq 0)]ds,$

which provides a Taylor-like expansion for non-smooth function. Then, we have

$\displaystyle\rho_{\tau}({Y_{i}}-{{X}}_{i}^{\top}\hat{\theta}_{n})$
$\displaystyle=\rho_{\tau}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta})+{{X}}_{i}^{\top}(\tilde{\theta}-\hat{\theta}_{n})[\tau-{\rm{1}}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta}\leq
0)]$ (45)
$\displaystyle+\int_{0}^{{{X}}_{i}^{\top}(\hat{\theta}_{n}-\tilde{\theta})}[{\rm{1}}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta}\leq
s)-{\rm{1}}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta}\leq 0)]ds$
$\displaystyle\geq\rho_{\tau}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta})-|{{X}}_{i}^{\top}(\tilde{\theta}-\hat{\theta}_{n})|\cdot|\tau-{\rm{1}}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta}\leq
0)|$
$\displaystyle-\left|\int_{0}^{{{X}}_{i}^{\top}(\hat{\theta}_{n}-\tilde{\theta})}[{\rm{1}}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta}\leq
s)-{\rm{1}}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta}\leq 0)]ds\right|$
$\displaystyle\geq\rho_{\tau}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta})-\max\\{\tau,1-\tau\\}|{{X}}_{i}^{\top}(\tilde{\theta}-\hat{\theta}_{n})|-|{{X}}_{i}^{\top}(\tilde{\theta}-\hat{\theta}_{n})|$
$\displaystyle[\text{By}~{}\eqref{eq:nets}]~{}$
$\displaystyle\geq\rho_{\tau}({Y_{i}}-{{X}}_{i}^{\top}\tilde{\theta})-\kappa\max\\{1+\tau,2-\tau\\}\|{{X}}_{i}\|_{{2}},~{}i=1,2,\cdots,n$
(46)

In below, let $l_{\tau}:=\max\\{1+\tau,2-\tau\\}$. The non-decreasing of
$\psi(\cdot)$ and (46) show that

$\displaystyle{{\hat{R}}_{\psi,{\rho_{\tau}},\alpha}}(\hat{\theta}_{n})=\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\hat{\theta}_{n})]\geq\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})-\alpha\kappa
l_{\tau}{\left\|{{X_{i}}}\right\|_{{2}}}].$ (47)

Next, we need to obtain the lower bound of the last term in (47) similar to
the proof of Lemma 15, but it includes covering number techniques. By the
lower bound in (3),

$\displaystyle~{}~{}~{}~{}{\rm{E}}{e^{-\sum_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})-\alpha\kappa
l_{\tau}{{\left\|{{X_{i}}}\right\|}_{{2}}}]}}$
$\displaystyle\leq{\rm{E\\{}}\prod\limits_{i=1}^{n}{[1-\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})+\alpha\kappa
l_{\tau}{{\left\|{{X_{i}}}\right\|}_{{2}}}+\frac{{{\alpha^{2}}}}{2}{[{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})-\kappa
l_{\tau}{{\left\|{{X_{i}}}\right\|}_{{2}}}]^{2}}]}{\rm{\\}}}$
$\displaystyle\leq\prod\limits_{i=1}^{n}{\rm{E}}[1-\alpha{\rho_{\tau}}(Y-{X^{\top}}\tilde{\theta})+\alpha\kappa
l_{\tau}{\left\|{{X_{i}}}\right\|_{{2}}}+\frac{{{\alpha^{2}}}}{2}{[{\rho_{\tau}}(Y-{X^{\top}}\tilde{\theta})-\kappa
l_{\tau}{\left\|{{X_{i}}}\right\|_{{2}}}]^{2}}]$ $\displaystyle\leq
e^{-\alpha\sum_{i=1}^{n}{\rm{E}}{\rho_{\tau}}(Y_{i}-{X_{i}^{\top}}\tilde{\theta})+\alpha\kappa
l_{\tau}\sum_{i=1}^{n}{\rm{E}}{\left\|{{X_{i}}}\right\|_{{2}}}+{\alpha^{2}}\sum_{i=1}^{n}\rho_{\tau}^{2}(Y_{i}-{X_{i}^{\top}}\tilde{\theta})+{(\alpha\kappa
l_{\tau})^{2}}\sum_{i=1}^{n}{\rm{E}}\left\|{{X_{i}}}\right\|_{{2}}^{2}}$
$\displaystyle\leq{e^{n\alpha\\{-{{R}_{{\rho_{\tau}}}^{n}}(\tilde{\theta})+\kappa
l_{\tau}{\rm{E}}{{\left\|{{X}}\right\|}_{{2}}}+\alpha{R_{\rho_{\tau}^{2}}^{n}}(\tilde{\theta})+\alpha{(\kappa
l_{\tau})^{2}}{\rm{E}}\left\|{{X}}\right\|_{{2}}^{2}\\}}},$ (48)

where the second last inequality stems from $1+x\leq e^{x}$ and $(a+b)^{2}\leq
2a^{2}+2b^{2}$.

By Markov’s inequality and the moment bound (A.7), we have for a fixed
$\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)$,

$\displaystyle
P\left\\{{\frac{{-1}}{{n\alpha}}\sum\limits_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})-\alpha\kappa
l_{\tau}{{\left\|{{X_{i}}}\right\|}_{{2}}}]\geq}\right.$
$\displaystyle\left.{\\{-{{R}_{{\rho_{\tau}}}^{n}}(\tilde{\theta})+\kappa
l_{\tau}{\rm{E}}{{\left\|{{X}}\right\|}_{2}}+\alpha{R_{\rho_{\tau}^{2}}^{n}}(\tilde{\theta})+\alpha{{(\kappa
l_{\tau})}^{2}}{\rm{E}}\left\|{{X}}\right\|_{{2}}^{2}\\}+{\frac{{\log(1/s)}}{{n\alpha}}}}\right\\}$
$\displaystyle\leq\frac{{{\rm{E}}{e^{-\sum_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})-\alpha\kappa
l_{\tau}{{\left\|{{X_{i}}}\right\|}_{{2}}}]}}}}{{{e^{n\alpha\\{-{{R}_{{\rho_{\tau}}}^{n}}(\tilde{\theta})+l_{\tau}{\rm{E}}{{\left\|{{X}}\right\|}_{{2}}}+\alpha{R_{\rho_{\tau}^{2}}^{n}}(\tilde{\theta})+\alpha{{(\kappa
l_{\tau})}^{2}}{\rm{E}}\left\|{{X}}\right\|_{{2}}^{2}\\}+\log(1/s)}}}}\leq
s\in(0,1).$ (49)

Note that the set $\mathcal{N}({\Theta},\kappa)$ has $N({\Theta},\kappa)$
elements. For all $\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)$, we consider
following union bound based on the single bound (A.7),
$s=\delta/N({\Theta},\kappa)$,

$\displaystyle
P\left[{\bigcup\limits_{\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)}{\left\\{{\frac{{-1}}{{n\alpha}}\sum\limits_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})-\alpha\kappa
l_{\tau}{{\left\|{{X_{i}}}\right\|}_{2}}]\geq}\right.}}\right.$
$\displaystyle\\{-{{R}_{{\rho_{\tau}}}^{n}}(\tilde{\theta})+\kappa
l_{\tau}{\rm{E}}{\left\|{{X}}\right\|_{2}}+\alpha{R_{\rho_{\tau}^{2}}^{n}}(\tilde{\theta})+\alpha{(\kappa
l_{\tau})^{2}}{\rm{E}}\left\|{{X}}\right\|_{2}^{2}\\}+\left.{\left.{\frac{{\log(1/s)}}{{n\alpha}}}\right\\}}\right]$
$\displaystyle\leq N({\Theta},\kappa)\cdot
P\left\\{{\frac{{-1}}{{n\alpha}}\sum\limits_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})-\alpha\kappa
l_{\tau}{{\left\|{{X_{i}}}\right\|}_{2}}]\geq}\right.$
$\displaystyle\\{-{{R}_{{\rho_{\tau}}}^{n}}(\tilde{\theta})+\kappa
l_{\tau}{\rm{E}}{\left\|{{X}}\right\|_{2}}+\alpha{R_{\rho_{\tau}^{2}}^{n}}(\tilde{\theta})+\alpha{(\kappa
l_{\tau})^{2}}{\rm{E}}\left\|{{X}}\right\|_{2}^{2}\\}+\left.{\frac{{\log(1/s)}}{{n\alpha}}}\right\\}\leq
N({\Theta},\kappa)s=:\delta.$ (50)

Continuing to bound (47),
$\forall~{}\tilde{\theta}\in\mathcal{N}({\Theta},\kappa)$, we have with
probability at least $1-\delta$ via (A.7)

$\displaystyle~{}~{}~{}~{}{{\hat{R}}_{\psi,{\rho_{\tau}},\alpha}}(\hat{\theta}_{n})=\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\hat{\theta}_{n})]\geq\frac{1}{{n\alpha}}\sum\limits_{i=1}^{n}\psi[\alpha{\rho_{\tau}}({Y_{i}}-X_{i}^{\top}\tilde{\theta})-\alpha\kappa
l_{\tau}{\left\|{{X_{i}}}\right\|_{2}}]$
$\displaystyle\geq{{R}_{{\rho_{\tau}}}^{n}}(\tilde{\theta})-[\kappa
l_{\tau}{\rm{E}}{\left\|{{X}}\right\|_{2}}+\alpha{R_{\rho_{\tau}^{2}}^{n}}(\tilde{\theta})+\alpha{(\kappa
l_{\tau})^{2}}{\rm{E}}\left\|{{X}}\right\|_{2}^{2}]+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta}]$
$\displaystyle\geq{{R}_{{\rho_{\tau}}}^{n}}(\tilde{\theta})-[\kappa
l_{\tau}{\rm{E}}{\left\|{{X}}\right\|_{2}}+\alpha\mathop{\sup}\limits_{\theta\in{\Theta}}{R_{\rho_{\tau}^{2}}^{n}}(\theta)+\alpha{(\kappa
l_{\tau})^{2}}{\rm{E}}\left\|{{X}}\right\|_{2}^{2}]+\frac{1}{{n\alpha}}\log\frac{N({\Theta},\kappa)}{\delta}].$
(51)

Let ${{\rm{E}}_{n}}f(X):=\frac{1}{n}\sum_{i=1}^{n}{{\rm{E}}f({X_{i}})}$ be the
average moment for the independent random variables $\\{{X_{i}}\\}_{i=1}^{n}$
in some space $\mathcal{X}$. Recall that
$\\{{Y_{i}}-X_{i}^{\top}\theta_{\tau}^{*}=:{\epsilon_{i}}\\}_{i=1}^{n}$ are
independent with densities ${f_{\varepsilon_{i}}}(t)$, where
${f_{\epsilon_{i}}}(t)$ is the conditional density function for $\epsilon_{i}$
given $X_{i}$. The quantile function gives the relation
${Q_{{Y_{i}|X_{i}}}}(\tau)=X_{i}^{\top}\theta_{n}^{*}\Leftrightarrow
P({Y_{i}}\leq X_{i}^{\top}\theta^{*}|X_{i})=\tau.$ To get the lower bound for
${{R}_{{\rho_{\tau}}}^{n}}(\tilde{\theta})$, we use Knight’s identity (45),

$\displaystyle~{}~{}~{}~{}{R_{{\rho_{\tau}}}^{n}}(\hat{\theta}_{n})={\rm{E}}_{n}[{\rho_{\tau}}(Y-{X^{\top}}\hat{\theta}_{n})-{\rho_{\tau}}(Y-{X^{\top}}\tilde{\theta})]+{R_{{\rho_{\tau}}}^{n}}(\tilde{\theta})$
$\displaystyle={R_{{\rho_{\tau}}}^{n}}(\tilde{\theta})+{\rm{E}}_{n}\\{X^{\top}(\hat{\theta}_{n}-\tilde{\theta})[\tau-{\rm{1}}({Y}-X^{\top}\tilde{\theta}\leq
0)]\\}$
$\displaystyle+{\rm{E}}_{n}\left({\mathop{\smallint}\nolimits_{0}^{X^{\top}(\hat{\theta}_{n}-\tilde{\theta})}[{\rm{1}}({Y}-X^{\top}\tilde{\theta}\leq
s)-{\rm{1}}({Y}-X^{\top}\tilde{\theta}\leq 0)]ds}\right)$
$\displaystyle\leq{R_{{\rho_{\tau}}}^{n}}(\tilde{\theta})+\max\\{\tau,1-\tau\\}{\rm{E}}|X^{\top}(\hat{\theta}_{n}-\tilde{\theta})|+{\rm{E}}|X^{\top}(\hat{\theta}_{n}-\tilde{\theta})|\leq{R_{{\rho_{\tau}}}^{n}}(\tilde{\theta})+l_{\tau}\kappa{\rm{E}}{\left\|X\right\|_{2}},$
(52)

i.e.
${R_{{\rho_{\tau}}}^{n}}(\tilde{\theta})\geq{R_{{\rho_{\tau}}}^{n}}(\hat{\theta}_{n})-l_{\tau}\kappa{\rm{E}}{\left\|X\right\|_{2}}$,
which implies (A.7) has a lower bound:

${{\hat{R}}_{\psi,{\rho_{\tau}},\alpha}}(\hat{\theta}_{n})\geq{{R}_{{\rho_{\tau}}}^{n}}(\hat{\theta}_{n})-[2\kappa
l_{\tau}{\rm{E}}{\left\|{{X}}\right\|_{2}}+\alpha\mathop{\sup}\limits_{\theta\in{\Theta}}{R_{\rho_{\tau}^{2}}^{n}}(\theta)+\alpha{(\kappa
l_{\tau})^{2}}{\rm{E}}\left\|{{X}}\right\|_{2}^{2}+\frac{1}{{n\alpha}}\log(\frac{N({\Theta},\kappa)}{\delta})]$