# The Dirichlet problem for prescribed curvature equations of $p$-convex
hypersurfaces

Weisong Dong  School of Mathematics, Tianjin University, 135 Yaguan Road,
Jinnan, Tianjin, China, 300354<EMAIL_ADDRESS>
###### Abstract.

In this paper, we study the Dirichlet problem for $p$-convex hypersurfaces
with prescribed curvature. We prove that there exists a graphic hypersurface
satisfying the prescribed curvature equation with homogeneous boundary
condition. An interior curvature estimate is also obtained.

_Keywords:_ Dirichlet problem; $C^{2}$ estimates; interior curvature
estimates, $p$-convex hypersurfaces.

##  1\. Introduction

In this paper, we are interested in finding a graphic hypersurface of
prescribed curvature on a bounded domain $\Omega\subset\mathbb{R}^{n}$.
Suppose that the hypersurface
$M_{u}:=\\{(x,u(x)):x\in\Omega\\}\subset\mathbb{R}^{n+1}$ is given by the
graph of a smooth function $u:\Omega\rightarrow\mathbb{R}$. Denote by
$\kappa[M_{u}]=(\kappa_{1},\cdots,\kappa_{n})$ the principal curvatures of
$M_{u}$ with respect to the downward unit normal of $M_{u}$. Given an integer
$p$, where $1\leqslant p\leqslant n$, it is interesting in geometry analysis
to find a graphic hypersurface $M_{u}$ with its principal curvatures
satisfying the following equation

(1.1) $\Pi_{1\leqslant i_{1}<\cdots<i_{p}\leqslant
n}(\kappa_{i_{1}}+\cdots+\kappa_{i_{p}})=f(x,u,\nu)\;\mbox{in}\;\Omega,$

and with the prescribed homogeneous boundary condition as below

(1.2) $u=0\;\mbox{on}\;\partial\Omega,$

where $\nu$ is the upward unit normal vector to the graphic hypersurface at
$X:=(x,u(x))$ and $f(x,z,\nu)>0$ is a smooth function defined on
$\Omega\times\mathbb{R}\times\mathbb{R}^{n}$. For convenience, we denote that
$\mathcal{M}_{u}^{p}:=\Pi_{1\leqslant i_{1}<\cdots<i_{p}\leqslant
n}(\kappa_{i_{1}}+\cdots+\kappa_{i_{p}})$. We also simply write
$\mathcal{M}_{u}$ when there is no ambiguity.

The equation (1.1) is a fully nonlinear partial differential equation for
$p<n$. It is natural to consider the problem in the class of $p$-convex
hypersurfaces so that the equation is elliptic. Recall that a $C^{2}$ regular
hypersurface $M_{u}$ is called (strictly) $p$-convex if $\kappa[M_{u}](x)$
satisfies, at each point $x\in\Omega$, that

$\kappa_{i_{1}}+\cdots+\kappa_{i_{p}}\geqslant\;(>)\;0,\>\forall\;1\leqslant
i_{1}<\cdots<i_{p}\leqslant n.$

Accordingly, we say that the $C^{2}$ function $u$ is _admissible_. The notion
of $p$-convexity can be dated back to Wu [51]. Sha in [41, 42] studied
Riemannian manifolds with $p$-convex boundaries and proved a complete
characterization for such Riemannian manifolds. Recently, there were some new
discoveries in $p$-convex geometry made by Harvey and Lawson. We refer the
reader to [23, 24] for these results and also to a nice survey [25].

There is a vast literature on the existence of closed Weingarten hypersurfaces
of codimension one in the Euclidean space. The Gausssian curvature case
corresponding to $p=1$ in (1.1) was studied by Oliker [36]. The mean curvature
equation corresponding to $p=n$ in (1.1) was studied by Bakelman-Kantor [1]
and Treibergs-Wei [47]. Caffarelli-Nirenberg-Spruck in [5] was concerned with
prescribed curvature equations of very general type, and it was then
generalized to Riemannian manifolds by Gerhardt in [13]. The existence of a
starshaped $(n-1)$-convex hypersurface of prescribed curvature (1.1) was
proved by Chu-Jiao [8]. The case $p\geqslant\frac{n}{2}$ was explored in a
subsequent work by the author in [10]. In complex settings, when $p=n-1$, the
operator appears in the Gauduchon conjecture which was solved by Székelyhidi-
Tosatti-Weinkove [44]. Some previous work on this topic can be found in
Tosatti-Weinkove [45, 46] and Fu-Wang-Wu [11, 12]. A counterpart to (1.1) that
was studied extensively is the higher order mean curvature equation

(1.3) $\sigma_{k}(\kappa)=f(X,\nu(X)),\;\forall X\in M,$

where $1\leqslant k\leqslant n$ and $\sigma_{k}$ is the $k$-th elementary
symmetric function. In warped product manifolds, Chen-Li-Wang [7] established
the curvature estimate for convex hypersurfaces satisfying (1.3). Jin-Li [30]
and Li-Oliker [33] investigated the starshaped hypersurfaces in hyperbolic and
elliptic space respectively. For the study in Minkowski space, see Li [31],
Wang-Xiao [50] and Ren-Wang-Xiao [40]. We also refer the reader to Guan-Zhang
[22] for a class of curvature equations arising from convex geometry, and to
Guan-Lin-Ma [20] and Guan-Li-Li [19] for hypersurfaces with prescribed
curvature measures.

The Dirichlet problem for Weingarten hypersurfaces has also attracted
attention of many authors, due to its role in solving a Plateau type problem
for locally convex Weingarten hypersurfaces. For its own sake, there still are
many works, such as Caffarelli-Nirenberg-Spruck [6], Lions [35], Ivochkina
[26, 27], Ivochkina-Lin-Trudinger [28], Lin-Trudinger [34] and Trudinger [48].
We remark that in [28] the Dirichlet problem with non-homogeneous boundary
data was considered for smooth solutions and in [48] it was investigated in
the context of viscosity solutions. In this paper, we are focused on the
Dirichlet problem of homogeneous boundary data. We shall return to the non-
homogeneous boundary case in future work. Now, we recall the $p$-convex cones
introduced by Harvey and Lawson [25].

######  Definition 1.1.

Let $p\in\\{1,\cdots,n\\}$. The $\mathcal{P}_{p}$ cone is a subset in
$\mathbb{R}^{n}$ with element $(\lambda_{1},\cdots,\lambda_{n})$ such that
$\lambda_{i_{1}}+\cdots+\lambda_{i_{p}}>0$ for all $1\leqslant
i_{1}<\cdots<i_{p}\leqslant n$. The $P_{p}$ cone of symmetric $n\times n$
matrices is defined as: $A\in P_{p}$ if the $n$-tuple of the eigenvalues of
matrix $A$ is in $\mathcal{P}_{p}$. We call $A$ is $p$-positive if $A\in
P_{p}$.

Denote by $(\kappa_{1}^{b},\cdots,\kappa_{n-1}^{b})$ the principal curvatures
of the boundary $\partial\Omega$ taken with respect to the exterior unit
normal to $\partial\Omega$ (see Section 14.6 in [14]). We call $\Omega$ is
_strictly convex_ if that
$(\kappa_{1}^{b}(x),\cdots,\kappa_{n-1}^{b}(x))\in\mathcal{P}_{1}\subset\mathbb{R}^{n-1}$
for every $x\in\partial\Omega$. Now, we state our main result.

######  Theorem 1.2.

Let $\Omega$ be a strictly convex bounded domain in $\mathbb{R}^{n}$ with
smooth boundary $\partial\Omega$ and $f=f(x,z,\nu)\in
C^{\infty}(\bar{\Omega}\times\mathbb{R}\times\mathbb{S}^{n})$ be a positive
function with $f_{z}\geqslant 0$. Suppose that $p\geqslant\frac{n}{2}$. Assume
that there is an admissible subsolution $\underline{u}\in C^{2}(\bar{\Omega})$
satisfying

(1.4) $\mathcal{M}_{\underline{u}}^{p}\geqslant
f(x,\underline{u},\underline{\nu})\;\mbox{in}\;\Omega\;\;\mbox{and}\;\;\underline{u}=0\;\mbox{on}\;\partial\Omega,$

where $\underline{\nu}$ denotes the upward unit normal vector to the graphic
hypersurface $M_{\underline{u}}$ at $\underline{X}:=(x,\underline{u}(x))$.
Then there exists a unique admissible solution $u\in C^{\infty}(\bar{\Omega})$
to (1.1) and (1.2).

######  Remark 1.3.

The endpoint case $p=n$ is excluded from our discussion, since this case
belongs to the realm of quasilinear elliptic equations. The case $p=1$ is also
well known. Jiao-Sun [29] proved the theorem for $p=n-1$ recently. It remains
an interesting question to prove the theorem for $1<p<\frac{n}{2}$.

The approach to solve the fully nonlinear elliptic equation (1.1) is the well-
known continuity method and degree theory. So the main task of this paper is
to establish _a priori_ estimates for admissible solutions up to second order
derivatives as higher order estimates follow from the Evans-Krylov theorem and
Schauder theory. By the assumption of the existence of a subsolution, $C^{0}$
and $C^{1}$ estimates are easy to derive. The most difficult part is to
establish _a priori_ $C^{2}$ estimates, including the global $C^{2}$ estimate
and the boundary $C^{2}$ estimate. Moreover, since $f$ depends on $\nu$, the
bad term $-Cv_{11}$ appears in the inequality (3.9) when applying the maximum
principle to the test function $G$ in Section 3, which arise from
differentiating the equation (1.1) twice (see (3.8)). This also cause trouble
for us to deal with third order terms. To overcome these difficulties, one can
impose various assumptions on $f$ as in [26, 27, 16, 17]. Due to the special
nature of the operator of equation (1.1), we use idea from [8, 10] to control
$-Cv_{11}$ by good terms without imposing any extra conditions on $f$. Then
the bad third order terms can be eliminated by the good term $-F^{1i,i1}$. See
Lemma 2.4 and Lemma 2.5. As to the boundary $C^{2}$ estimate, the issue is to
estimate the tangential-normal derivatives of solutions. Fortunately, we can
use idea from Ivochkina [26] to tackle this. Actually, we have enough good
third order terms to derive an interior curvature estimate.

######  Theorem 1.4.

Let $\Omega$ be a bounded domain in $\mathbb{R}^{n}$. Suppose that $u\in
C^{4}(\Omega)\cap C^{2}(\bar{\Omega})$ is an admissible solution to (1.1) and
(1.2). Then there is a positive constant $C$ and $\beta$ such that the second
fundamental form $h$ of $M_{u}$ satisfies

(1.5) $\sup_{\Omega}(-u)^{\beta}|h|\leqslant C,$

where $C$ and $\beta$ depend on $n$, $p$, $|u|_{C^{1}}$, $\inf\tilde{f}$ and
$|f|_{C^{2}}$.

Some applications of the interior curvature estimate can be found in Sheng-
Urbas-Wang [43]. The above theorem holds for $p\geqslant\frac{n}{2}$ if
$f\equiv f(x,u,\nu)$ and for $1\leqslant p\leqslant n$ if $f\equiv f(x,u)$ as
explained below. The assumption $p\geqslant\frac{n}{2}$ in Theorem 1.2 is only
used to prove Lemma 2.4 to control $-Cv_{11}$ in the inequality (3.10) (and in
(3.21) for the interior curvature estimate). If the right hand side function
$f$ in the equation (1.1) does not depend on $\nu$, the term $-Cv_{11}$ does
not appear in (3.10) (also not appear in (3.21)), for which we refer the
reader to examine the inequalities (3.8) and (3.9). Therefore, we have the
following result.

######  Theorem 1.5.

Let $\Omega$ be a strictly convex bounded domain in $\mathbb{R}^{n}$ with
smooth boundary $\partial\Omega$ and $f=f(x,z)\in
C^{\infty}(\bar{\Omega}\times\mathbb{R})$ be a positive function with
$f_{z}\geqslant 0$. Suppose $1\leqslant p\leqslant n$. Assume that there is an
admissible subsolution $\underline{u}\in C^{2}(\bar{\Omega})$ satisfying

(1.6) $\mathcal{M}_{\underline{u}}^{p}\geqslant
f(x,\underline{u})\;\mbox{in}\;\Omega\;\;\mbox{and}\;\;\underline{u}=0\;\mbox{on}\;\partial\Omega,$

where $\underline{\nu}$ denotes the upward unit normal vector to the graphic
hypersurface $M_{\underline{u}}$ at $\underline{X}:=(x,\underline{u}(x))$.
Then there exists a unique admissible solution $u\in C^{\infty}(\bar{\Omega})$
to (1.1) and (1.2).

Before we end the introduction, let’s review a conjecture on curvature
estimates for hypersurfaces with prescribed higher order mean curvature as it
is closely related to our topic. The Gaussian curvature and the mean curvature
correspond to $k=n$ and $k=1$ in (1.3) respectively. Ren-Wang conjectured that
the curvature estimate holds for $k>\frac{n}{2}$ in [38, 39]. While the case
$k=n$ proved a long time ago by Caffarelli-Nirenberg-Spruck [3], it was only
recently that Ren-Wang [37, 38] proved the conjecture for $k\geqslant n-2$.
Guan-Ren-Wang in [21] established the curvature estimate when $k=2$, and in
the same paper they also proved the curvature estimate for convex solutions of
equation (1.3) when $3\leqslant k\leqslant n$.

The paper is organized as follows. In Section 2, we recall some properties of
the operator from [9] and prove some key lemmas which are crucial to the
global and boundary $C^{2}$ estimate. Then, we prove the global and boundary
$C^{2}$ estimate in Section 3 and Section 4 respectively. In Section 5, we
provide a proof of the gradient estimate.

Acknowledgements: The author is supported by the National Natural Science
Foundation of China, No. 11801405\.

##  2\. Preliminaries

Let $M$ be a hypersurface in $\mathbb{R}^{n+1}$ given by the graph of
$u:\Omega\rightarrow\mathbb{R}$. Then the induced metric $g$ and second
fundamental form $h$ of $M$ are given by

$g_{ij}=\delta_{ij}+u_{i}u_{j}\;\mbox{and}\;h_{ij}=\frac{u_{ij}}{\sqrt{1+|Du|^{2}}}.$

The upward unit normal vector is

$\nu=\frac{(-Du,1)}{\sqrt{1+|Du|^{2}}}.$

The principal curvatures $\kappa=(\kappa_{1},\cdots,\kappa_{n})$ of the
hypersurface $M$ are the eigenvalues of $h_{ij}$ with respect to $g_{ij}$,
i.e. the eigenvalues of the following matrix

(2.1) $\frac{1}{w}\Big{(}I-\frac{Du\otimes Du}{w^{2}}\Big{)}D^{2}u,$

where $w=\sqrt{1+|Du|^{2}}$ and $g^{ij}=\delta_{ij}-\frac{u_{i}u_{j}}{w^{2}}$
is the inverse of $g_{ij}$. Note that the above matrix generally is not
symmetric. As in [6], we consider the following matrix

(2.2) $a_{ij}=\frac{1}{w}\gamma^{ik}u_{kl}\gamma^{lj},$

where $\gamma^{ik}=\delta_{ik}-\frac{u_{i}u_{k}}{w(1+w)}$ is the inverse
matrix of $\gamma_{ij}=\delta_{ij}+\frac{u_{i}u_{j}}{1+w}$ and $\gamma_{ij}$
is the square root of $g_{ij}$.

For convenience, we introduce the following notations

$F(a_{ij}):=F(\kappa)=\Pi_{1\leqslant i_{1}<\cdots<i_{p}\leqslant
n}(\kappa_{i_{1}}+\cdots+\kappa_{i_{p}})\;\mbox{and}\;\tilde{F}=F^{\frac{1}{C_{n}^{p}}},$

where $C_{n}^{p}=\frac{n!}{p!(n-p)!}$. Equation (1.1) then can be written as

(2.3) $\tilde{F}(a_{ij}):=\tilde{F}(\kappa)=\tilde{f}(x,u,\nu),$

where $\kappa=(\kappa_{1},\cdots,\kappa_{n})$ and
$\tilde{f}=f^{\frac{1}{C_{n}^{p}}}$. It is also convenient to rewrite the
equation as below

(2.4) $G(D^{2}u,Du):=\tilde{F}(a_{ij})=\tilde{f},$

where $a_{ij}$ is defined in (2.2). Denote

$F^{ij}=\frac{\partial F}{\partial
a_{ij}},\;F^{ij,kl}=\frac{\partial^{2}F}{\partial a_{ij}\partial
a_{kl}},\;G^{ij}=\frac{\partial G}{\partial u_{ij}},\;G^{i}=\frac{\partial
G}{\partial u_{i}}.$

We remark that $G$ is concave and homogeneous one with respect to $D^{2}u$
since $\tilde{F}$ is concave and homogeneous one with respect to $a_{ij}$. See
[9]. The equation is elliptic as the matrix
$\\{\frac{\partial\tilde{F}}{\partial a_{ij}}\\}$ is positive definite for
$\\{a_{ij}\\}\in P_{p}$. It is easy to verify that

$\frac{1}{w}\sum\tilde{F}^{ii}\geqslant\sum
G^{ii}\geqslant\frac{1}{w^{3}}\sum\tilde{F}^{ii}.$

By a straightforward calculation, one get the following lemma. See formula
(2.21) in Lemma 2.3 of [18].

######  Lemma 2.1 ([18]).

We have

(2.5)
$G^{s}=-\frac{u_{s}}{w^{2}}\sum_{i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}-\frac{2}{w(1+w)}\sum_{ij}\frac{\partial\tilde{F}}{\partial
a_{ij}}a_{it}(wu_{t}\gamma^{sj}+u_{j}\gamma^{ts}).$

Next, we calculate the derivatives of $F$. The calculations are carried out at
a point $X_{0}$ on the hypersurface $M_{u}$ where $a_{ij}$ is diagonal. Note
that $F^{ij}$ is also diagonal at $X_{0}$ and we have the following formulas

$F^{kk}=\frac{\partial
F}{\partial\kappa_{k}}=\sum_{k\in\\{i_{1},\cdots,i_{p}\\}}\frac{F(\kappa)}{\kappa_{i_{1}}+\cdots+\kappa_{i_{p}}},$

for which we refer to Lemma 1.10 in [9]. We have formulas for the second order
derivatives of $F$ at $X_{0}$ as below

$F^{kk,ll}=\frac{\partial^{2}F}{\partial\kappa_{k}\partial\kappa_{l}}=\sum_{\begin{subarray}{c}k\in\\{i_{1},\cdots,i_{p}\\}\\\
l\in\\{j_{1},\cdots,j_{p}\\}\\\
\\{i_{1},\cdots,i_{p}\\}\neq\\{j_{1},\cdots,j_{p}\\}\end{subarray}}\frac{F(\kappa)}{(\kappa_{i_{1}}+\cdots+\kappa_{i_{p}})(\kappa_{j_{1}}+\cdots+\kappa_{j_{p}})},$

and, for $k\neq r$,

$F^{kr,rk}=\frac{F^{kk}-F^{rr}}{\kappa_{k}-\kappa_{r}}=-\sum_{\begin{subarray}{c}k\notin\\{i_{1},\cdots,i_{p}\\}\ni
r\\\ r\notin\\{j_{1},\cdots,j_{p}\\}\ni k\\\
\\{i_{1},\cdots,i_{p}\\}\backslash\\{r\\}=\\{j_{1},\cdots,j_{p}\\}\backslash\\{k\\}\end{subarray}}\frac{F(\kappa)}{(\kappa_{i_{1}}+\cdots+\kappa_{i_{p}})(\kappa_{j_{1}}+\cdots+\kappa_{j_{p}})}.$

Otherwise, we have $F^{ij,kl}=0$. See Lemma 1.12 in [9] for the above
formulas, or see [10]. These formulas can also be easily obtained from Theorem
5.5 in [2]. The following properties of the function $F$ which are very
similar to the properties of $\sigma_{k}$ were proved by Dinew [9].

######  Lemma 2.2 ([9]).

Suppose that the diagonal matrix $A=diag(\lambda_{1},\cdots,\lambda_{n})$
belongs to $P_{p}$ and that $\lambda_{1}\geqslant\cdots\geqslant\lambda_{n}$.
Then,

  * (1)

$\tilde{F}^{11}(A)\lambda_{1}\geqslant\frac{1}{n}\tilde{F}(A);$

  * (2)

$\sum_{k=1}^{n}\tilde{F}^{kk}(A)\geqslant p;$

  * (3)

$\sum_{k=1}^{n}F^{kk}(A)\lambda_{k}=C_{n}^{p}F(A);$

  * (4)

there is a constant $\theta=\theta(n,p)$ such that
$F^{jj}(A)\geqslant\theta\sum F^{ii}$ for all $j\geqslant n-p+1.$

The proof of the above Lemma can also be found in the Appendix of [10]. Next,
we prove a few lemmas which were used to derive second order estimates.

######  Lemma 2.3.

For every $C>0$ and every compact set $K\subset\mathcal{P}_{p}$, there exists
a number $R=R(C,K)$ such that the inequality

(2.6) $F(\lambda_{1},\cdots,\lambda_{n}+R)\geqslant C$

holds for all $\lambda=(\lambda_{1},\cdots,\lambda_{n})\in K$.

###### Proof.

Since $K$ is compact, there is a positive constant $\varrho=C(K)$ such that,
for all $\lambda\in K$ and
$\\{i_{1},\cdots,i_{p}\\}\subset\\{1,\cdots,n-1\\}$,

$\lambda_{i_{1}}+\cdots+\lambda_{i_{p}}\geqslant\varrho.$

Therefore, for all $\lambda\in K$, we have

(2.7) $\displaystyle\;F(\lambda_{1},\cdots,\lambda_{n}+R)$ $\displaystyle=$
$\displaystyle\;\Pi_{1\leqslant i_{1}<\cdots<i_{p-1}\leqslant
n-1}(\lambda_{i_{1}}+\cdots+\lambda_{i_{p-1}}+\lambda_{n}+R)$
$\displaystyle\;\times\Pi_{1\leqslant i_{1}<\cdots<i_{p}\leqslant
n-1}(\lambda_{i_{1}}+\cdots+\lambda_{i_{p}})$ $\displaystyle\geqslant$
$\displaystyle\;R^{C_{n-1}^{p-1}}\times\Pi_{1\leqslant
i_{1}<\cdots<i_{p}\leqslant n-1}(\lambda_{i_{1}}+\cdots+\lambda_{i_{p}})$
$\displaystyle\geqslant$
$\displaystyle\;R^{C_{n-1}^{p-1}}\varrho^{C_{n-1}^{p}},$

which goes to infinity when $R$ goes to infinity. Thus, the lemma is proved. ∎

The following two lemmas were essentially proved in [10] (see Lemma 3.2 and
Lemma 3.4), which were used to deal with the bad third order terms. Here, we
restate them in the following form.

######  Lemma 2.4.

Suppose $p\geqslant\frac{n}{2}$ and $A=diag(\lambda_{1},\cdots,\lambda_{n})\in
P_{p}$ with $\lambda_{1}\geqslant\cdots\geqslant\lambda_{n}$. If
$\lambda_{n}\geqslant-\delta\lambda_{1}$ where
$0<\delta\leqslant\frac{1}{2(p-1)}$, we have

$C\sum F^{ii}(A)\geqslant\lambda_{1}$

for $C$ sufficiently large.

###### Proof.

We divide the proof into two cases.

Case 1. Suppose
$\lambda_{n-p+1}+\lambda_{n-p+2}+\cdots+\lambda_{n}<\frac{1}{\lambda_{1}}$.
Since

$F^{nn}\geqslant\frac{F(\kappa)}{\lambda_{n-p+1}+\lambda_{n-p+2}+\cdots+\lambda_{n}},$

we see that

$F^{nn}\geqslant c_{0}\lambda_{1}$

for some $c_{0}>0$ depending on $\inf f$. The lemma follows for sufficiently
large $C$ depending on $\inf f$.

Case 2. Suppose
$\lambda_{n-p+1}+\lambda_{n-p+2}+\cdots+\lambda_{n}\geqslant\frac{1}{\lambda_{1}}$.
For a fixed $(p-1)$-tuple $2\leqslant i_{1}<\cdots<i_{p-1}\leqslant n$, we
have

$\lambda_{1}+\lambda_{i_{1}}+\cdots+\lambda_{i_{p-1}}\geqslant\frac{\lambda_{1}}{2}$

by our assumption. Hence, we have

$\displaystyle F^{nn}\geqslant$
$\displaystyle\;\mathop{\Pi}\limits_{2\leqslant i_{1}<\cdots<i_{p-1}\leqslant
n}(\lambda_{1}+\lambda_{i_{1}}+\cdots+\lambda_{i_{p-1}})$
$\displaystyle\times\mathop{\Pi}\limits_{\begin{subarray}{c}2\leqslant
i_{1}<\cdots<i_{p}\leqslant n\\\
(i_{1},\cdots,i_{p})\neq(n-p+1,\cdots,n)\end{subarray}}(\lambda_{i_{1}}+\lambda_{i_{2}}+\cdots+\lambda_{i_{p}})$
$\displaystyle\geqslant$
$\displaystyle\;[\frac{\lambda_{1}}{2}]^{C_{n-1}^{p-1}}[\frac{1}{\lambda_{1}}]^{C_{n-1}^{p}-1}.$

For $p\geqslant\frac{n}{2}$, a direct calculation shows that

$C_{n-1}^{p-1}-C_{n-1}^{p}=\frac{(n-1)\cdots(n-p+1)}{(p-1)!}\Big{(}1-\frac{n-p}{p}\Big{)}\geqslant
0.$

Therefore, we obtain

$2^{C_{n-1}^{p-1}}F^{nn}\geqslant\lambda_{1}.$

∎

######  Lemma 2.5.

Suppose $A=diag(\lambda_{1},\cdots,\lambda_{n})\in P_{p}$ with
$\lambda_{1}\geqslant\cdots\geqslant\lambda_{n}$. For any given
$1>\epsilon>0$, there is a $\delta=\delta(\epsilon)>0$ such that, if
$\lambda_{n}\geqslant-\delta\lambda_{1}$, the inequality

$-2F^{1i,i1}(A)+(1+\epsilon)\frac{F^{11}(A)}{\lambda_{1}}\geqslant(1+\epsilon)\frac{F^{ii}(A)}{\lambda_{1}}$

holds for $i=2,3,\cdots n$.

###### Proof.

Observe that
$\frac{F^{ii}}{\lambda_{1}}=\frac{\lambda_{1}-\lambda_{i}}{\lambda_{1}}(-F^{1i,i1})+\frac{F^{11}}{\lambda_{1}}$,
from which we obtain that

$\frac{F^{ii}}{\lambda_{1}}\leqslant-F^{1i,i1}+\frac{F^{11}}{\lambda_{1}}\;\mbox{for}\;\lambda_{i}\geqslant
0.$

Since $-\lambda_{n}\leqslant\delta\lambda_{1}$, we have

$\frac{F^{ii}}{\lambda_{1}}\leqslant-(1+\delta)F^{1i,i1}+\frac{F^{11}}{\lambda_{1}}\;\mbox{for}\;\lambda_{i}\leqslant
0.$

By the above two inequalities, we get the desired inequality for sufficiently
small $\delta$. ∎

##  3\. Global and interior $C^{2}$ estimates

Let $\langle\cdot,\cdot\rangle$ be the inner product in $\mathbb{R}^{n+1}$ and
$\\{e_{1},\cdots,e_{n},e_{n+1}\\}$ be the standard basis of
$\mathbb{R}^{n+1}$. Then $X=(x,u(x))=x_{1}e_{1}+\cdots+x_{n}e_{n}+u(x)e_{n+1}$
and $|X|^{2}=\langle X,X\rangle=\sum_{i}x_{i}^{2}+u(x)^{2}$.

###  3.1. Global $C^{2}$ estimates

Now we prove the global $C^{2}$ estimate. The main idea which is from [10] is
to control the bad third order terms.

######  Theorem 3.1.

Let $\Omega$ be a bounded domain in $\mathbb{R}^{n}$. Suppose that $u\in
C^{4}(\Omega)\cap C^{2}(\bar{\Omega})$ is an admissible solution to (1.1).
Then there is a positive constant $C$ such that

(3.1) $\sup_{\Omega}|D^{2}u|\leqslant C(1+\sup_{\partial\Omega}|D^{2}u|),$

where $C$ depends on $n$, $p$, $|u|_{C^{1}}$, $\inf f$ and $|f|_{C^{2}}$.

###### Proof.

By the relation between $h_{ij}$ and $u_{ij}$, to prove the global $C^{2}$
estimate it is equivalent to prove the curvature estimate. Let $a$ be a
positive constant depending on $|Du|_{C^{0}}$ such that $\frac{1}{w}\geqslant
2a$. Consider the following test function

$G(x,\xi)=(\frac{1}{w}-a)^{-1}e^{\frac{A}{2}|X|^{2}}h_{\xi\xi},$

where $X\in M_{u}$, $\xi\in T_{X}M_{u}$ is a unit vector and $A$ is a large
constant to be determined later. Suppose the maximum of $G$ is achieved at a
point $X_{0}=(x_{0},u(x_{0}))\in M_{u}$ and $\xi_{0}\in T_{X_{0}}M_{u}$. We
choose new orthonormal vectors
$\\{\epsilon_{1},\cdots,\epsilon_{n},\epsilon_{n+1}\\}$ at $X_{0}$ such that
$\xi_{0}=\epsilon_{1}$ and $\nu(X_{0})=\epsilon_{n+1}$. Denote
$x_{0}=(x_{1}^{0},\cdots,x_{n}^{0})$. As in [6], we can represent the
hypersurface near $X_{0}$ by the tangential coordinates $y_{1},\cdots,y_{n}$
and $v(y)$:

$X=\sum_{j=1}^{n}x^{0}_{j}e_{j}+u(x_{0})e_{n+1}+\sum_{j=1}^{n}y_{j}\epsilon_{j}+v(y)\epsilon_{n+1}.$

Thus, we have $v(0)=0$ and $\nabla v(0)=0$. Using the same notations as in
[6], set $\omega=(1+|\nabla v|^{2})^{1/2}$. Then the principal curvature in
the $\epsilon_{1}$ direction is

$\kappa_{1}=\frac{v_{11}}{(1+v_{1}^{2})\omega}.$

Under the new coordinates, we denote the point $X\in M_{u}$ near $X_{0}$ by
$Y=\sum_{j=1}^{n}y_{j}\epsilon_{j}+v(y)\epsilon_{n+1}$. The normal in the new
coordinates still denoted by $\nu$ is given as

$\nu(Y)=-\frac{1}{\omega}\sum_{j=1}^{n}v_{j}\epsilon_{j}+\frac{1}{\omega}\epsilon_{n+1}.$

Recall that
$\epsilon_{n+1}=\frac{1}{w(x_{0})}(-u_{1}(x_{0}),\cdots,-u_{n}(x_{0}),1)$.
Then,

(3.2) $\frac{1}{w}=\nu(Y)\cdot e_{n+1}=\frac{1}{\omega
w(x_{0})}-\frac{1}{\omega}\sum_{j=1}^{n}a_{j}v_{j},$

where $a_{j}=\epsilon_{j}\cdot e_{n+1}$. It is easy to see that
$\sum_{j}a_{j}^{2}\leqslant 1$.

We compute that $|Y|^{2}=\sum_{j}y_{j}^{2}+v(y)^{2}$ and
$|X|^{2}=|X_{0}|^{2}+|Y|^{2}+2\langle X_{0},Y\rangle$. At the point $y=0$ the
function

(3.3)
$G(y)=(1/w-a)^{-1}e^{\frac{A}{2}|X|^{2}}\frac{v_{11}}{(1+v_{1}^{2})\omega}$

attains its maximum that is equal to $G(x_{0},\xi_{0})$. At this point, we
also have $v_{1j}=0$ for $j=2,\cdots,n$ since the $y_{1}$ axe lines along the
principal direction $\epsilon_{1}$. By a rotation of coordinates
$y_{2},\cdots,y_{n}$, we may assume that $v_{ij}(0)$ is diagonal. By (2.2) and
$\nabla v(0)=0$, we find that $a_{il}=v_{il}=\delta_{il}\kappa_{i}$ at the
origin. This implies that $\tilde{F}^{ij}$ is diagonal. Without loss of
generality, we assume that

$\kappa_{1}\geqslant\kappa_{2}\geqslant\cdots\geqslant\kappa_{n}.$

At the origin, by direct calculations (see (2.12) in [6]), we see that
$\frac{\partial a_{il}}{\partial y_{j}}=a_{ilj}=v_{ilj}$ and

$\frac{\partial^{2}a_{il}}{\partial
y_{1}^{2}}=a_{il11}=v_{il11}-v_{11}^{2}(v_{il}+\delta_{i1}v_{1l}+\delta_{1l}v_{1i}).$

Recall that $v(0)=0$ and $\nabla v(0)=0$. Differentiating $\log G(y)$ at $y=0$
twice yields that, for $1\leqslant i\leqslant n$,

(3.4) $0=-\frac{(1/w)_{i}}{1/w-a}+A(y_{i}+\langle
X_{0},\epsilon_{i}\rangle)+\frac{v_{11i}}{v_{11}}-\frac{2v_{1}v_{1i}}{1+v_{1}^{2}}-\frac{\omega_{i}}{\omega}$

and

(3.5)
$0\geqslant-\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}_{i}+A+\frac{v_{11ii}}{v_{11}}-\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-2v_{1i}^{2}-v_{ii}^{2}.$

Contracting (3.5) with $\tilde{F}^{ii}$, we get

(3.6) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{\tilde{F}^{ii}v_{11ii}}{v_{11}}-\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-\frac{\tilde{F}^{ii}(1/w)_{ii}}{1/w-a}$
$\displaystyle\;+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}-2\tilde{F}^{11}v_{11}^{2}-\tilde{F}^{ii}v_{ii}^{2}+A\sum\tilde{F}^{ii}.$

By (3.2), differentiate $1/w$ with respect to $y_{i}$ to obtain

$(1/w)_{i}=-a_{i}v_{ii}\;\mbox{and}\;(1/w)_{ii}=-\sum_{j}a_{j}v_{jii}-\frac{v_{ii}^{2}}{w(x_{0})}$

for all $1\leqslant i\leqslant n$. Combining the above formulas with the
derivatives of $a_{il}$, we therefore have from (3.6) that

(3.7) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{\tilde{F}^{ii}a_{ii11}}{v_{11}}+v_{11}\tilde{F}^{ii}v_{ii}-\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}+\frac{a_{j}\tilde{F}^{ii}a_{iij}}{1/w-a}$
$\displaystyle\;+\frac{\tilde{F}^{ii}v_{ii}^{2}}{1-aw}+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}-\tilde{F}^{ii}v_{ii}^{2}+A\sum\tilde{F}^{ii}.$

By differentiating the equation
$\tilde{F}(a_{ij})=\tilde{f}(X,\nu(X))=\tilde{f}(Y,\nu(Y))$ twice, we obtain

$\tilde{F}^{ii}a_{iij}=(\tilde{f})_{j}\geqslant-C-Cv_{11}$

and

(3.8) $\displaystyle\tilde{F}^{ii}a_{ii11}=$
$\displaystyle\;-\tilde{F}^{ij,kl}a_{ij1}a_{kl1}+(\tilde{f})_{11}$
$\displaystyle\geqslant$
$\displaystyle\;-\tilde{F}^{ij,kl}a_{ij1}a_{kl1}-Cv_{11}^{2}-\sum_{j}\tilde{f}_{\nu_{j}}\frac{v_{j11}}{\omega}$

for sufficiently large $v_{11}$, where $C$ depends on $|f|_{C^{2}}$ and
$|u|_{C^{1}}$. Substituting the above inequalities into (3.7), we arrive at

(3.9) $\displaystyle 0\geqslant$
$\displaystyle\;-\frac{2}{v_{11}}\sum_{i\geqslant
2}\tilde{F}^{1i,i1}a_{1i1}a_{i11}-Cv_{11}-\sum_{j}\tilde{f}_{\nu_{j}}\frac{v_{j11}}{v_{11}}-C$
$\displaystyle\;-\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\frac{aw}{1-aw}\tilde{F}^{ii}v_{ii}^{2}+A\sum\tilde{F}^{ii},$

provided that $v_{11}$ is large enough, where we used by Lemma 2.2 (3) that
$\tilde{F}^{ii}a_{ii}=\tilde{f}$. Using (3.4), we further derive that

(3.10) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{aw}{1-aw}\tilde{F}^{ii}v_{ii}^{2}-\frac{2}{v_{11}}\sum_{i\geqslant
2}\tilde{F}^{1i,i1}a_{1i1}a_{i11}-CA$
$\displaystyle\;-\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+A\sum\tilde{F}^{ii}-Cv_{11}$

for large $A$.

We now begin to deal with the third order terms. The following lemma ensures
us that in what follows we can assume that
$\kappa_{n}\geqslant-\delta\kappa_{1}$, where
$\delta\leqslant\frac{1}{2(p-1)}$ is a small positive constant to be
determined.

######  Lemma 3.2.

Suppose that $\kappa_{n}\leqslant-\delta\kappa_{1}$. Then, we can derive an
upper bound for $\kappa_{1}$ from (3.10).

###### Proof.

By the critical equation (3.4) and the Cauchy-Schwarz inequality, we have

$\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}\leqslant(1+\varepsilon)\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\Big{(}1+\frac{1}{\varepsilon}\Big{)}A^{2}\tilde{F}^{ii}(y_{i}+\langle
X_{0},\epsilon_{i}\rangle)^{2}$

for any $\varepsilon>0$. From (3.10) and $-\tilde{F}^{1i,i1}\geqslant 0$, we
see that

(3.11) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{aw}{1-aw}\tilde{F}^{ii}v_{ii}^{2}-\varepsilon\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}-\frac{CA^{2}}{\varepsilon}\sum\tilde{F}^{ii}-Cv_{11}$

for $v_{11}$ sufficiently large. Using $(1/w)_{i}=-a_{i}v_{ii}$ and choosing
$\varepsilon$ sufficiently small, we obtain from (3.11) that

$\displaystyle 0\geqslant$
$\displaystyle\;\frac{a}{2}\tilde{F}^{ii}v_{ii}^{2}-\frac{CA^{2}}{\varepsilon}\sum\tilde{F}^{ii}-Cv_{11}$
$\displaystyle\geqslant$
$\displaystyle\;\frac{a\theta}{2}\big{(}\sum\tilde{F}^{ii}\big{)}\delta^{2}\kappa_{1}^{2}-\frac{CA^{2}}{\varepsilon}\sum\tilde{F}^{ii}-Cv_{11},$

where we used Lemma 2.2 (4) in the second inequality. This implies an upper
bound for $\kappa_{1}$ since $\sum\tilde{F}^{ii}\geqslant p$.

∎

For an $\epsilon>0$ sufficiently small, we choose $\delta=\delta(\epsilon)$
such that Lemma 2.5 holds. Now we assume
$\kappa_{n}\geqslant-\delta\kappa_{1}$. By Lemma 2.4, the inequality (3.10)
becomes

(3.12) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{aw}{1-aw}\tilde{F}^{ii}v_{ii}^{2}-\frac{2}{v_{11}}\sum_{i\geqslant
2}\tilde{F}^{1i,i1}v_{1i1}v_{i11}-CA$
$\displaystyle\;-\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\frac{A}{2}\sum\tilde{F}^{ii}$

as long as $A$ is sufficiently large. For any given $\epsilon>0$, by Lemma
2.5, (3.12) becomes

(3.13) $\displaystyle 0\geqslant$
$\displaystyle\,a\tilde{F}^{ii}v_{ii}^{2}+\sum_{i\geqslant
2}\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-CA$
$\displaystyle\;-(1+\epsilon)\sum_{i\geqslant
2}\tilde{F}^{11}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\frac{A}{2}\sum\tilde{F}^{ii}.$

By the critical equation (3.4) and the Cauchy-Schwarz inequality, we see that

(3.14)
$\Big{(}\frac{v_{111}}{v_{11}}\Big{)}^{2}\leqslant(1+\epsilon)\Big{(}\frac{(1/w)_{1}}{1/w-a}\Big{)}^{2}+\frac{CA^{2}}{\epsilon}$

and, for $i\geqslant 2$,

(3.15)
$(1+\epsilon)\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}\leqslant(1+\epsilon)^{2}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\frac{CA^{2}}{\epsilon}.$

Substituting (3.14) and (3.15) into (3.13), we have

(3.16) $\displaystyle 0\geqslant$
$\displaystyle\,a\tilde{F}^{ii}v_{ii}^{2}-(1+\epsilon)\tilde{F}^{11}\Big{(}\frac{(1/w)_{1}}{1/w-a}\Big{)}^{2}-\frac{CA^{2}}{\epsilon}\tilde{F}^{11}-CA$
$\displaystyle\;-(1+3\epsilon)\tilde{F}^{11}\sum_{i\geqslant
2}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\frac{A}{2}\sum\tilde{F}^{ii}.$

Observe that $|v_{ii}|\leqslant nv_{11}$ and
$\tilde{F}^{11}\leqslant\cdots\leqslant\tilde{F}^{nn}$. We further derive that

(3.17) $\displaystyle 0\geqslant$
$\displaystyle\,a\tilde{F}^{ii}v_{ii}^{2}-\frac{3n^{2}\epsilon}{a^{2}}\tilde{F}^{11}v_{11}^{2}-\frac{CA^{2}}{\epsilon}\tilde{F}^{11}-CA.$

Choose $\epsilon\leqslant\frac{a^{3}}{6n^{2}}$. We then derive an upper bound
for $v_{11}$ by Lemma 2.2 (1).

∎

###  3.2. Interior $C^{2}$ estimates

As in Sheng-Urbas-Wang [43], we actually can derive an interior curvature
estimate for equation (1.1) with $p\geqslant\frac{n}{2}$ by the same argument
as that of Theorem 1.3 in [10]. Furthermore, if $f$ does not depend on $\nu$,
the interior estimate holds for all $p$.

######  Theorem 3.3.

Let $\Omega$ be a bounded domain in $\mathbb{R}^{n}$. Suppose that $u\in
C^{4}(\Omega)\cap C^{2}(\bar{\Omega})$ is an admissible solution to (1.1) and
(1.2). Then there is a positive constant $C$ and $\beta$ such that

(3.18) $\sup_{\Omega}(-u)^{\beta}|h|\leqslant C,$

where $C$ and $\beta$ depend on $n$, $p$, $|u|_{C^{1}}$, $\inf\tilde{f}$ and
$|f|_{C^{2}}$.

###### Proof.

We consider the test function
$G(x,\xi)=(-u)^{\beta}(\frac{1}{w}-a)^{-1}e^{\frac{A}{2}|X|^{2}}h_{\xi\xi}$.
Assume its maximum is attained at a point $X_{0}=(x_{0},u(x_{0}))\in M_{u}$
and $\xi_{0}\in T_{X_{0}}M_{u}$. As before, after choosing new coordinates at
$X_{0}$, the maximum of

$G(y)=(-u)^{\beta}(\frac{1}{w}-a)^{-1}e^{\frac{A}{2}|X|^{2}}\frac{v_{11}}{(1+v_{1}^{2})\omega}$

is achieved at $y=0$ and $a_{il}=v_{il}=\kappa_{i}\delta_{il}$. Without loss
of generality, we also assume that
$\kappa_{1}\geqslant\kappa_{2}\geqslant\cdots\geqslant\kappa_{n}$; then
$\tilde{F}^{11}\leqslant\tilde{F}^{22}\leqslant\cdots\leqslant\tilde{F}^{nn}$.
Note that $u=X\cdot e_{n+1}=X_{0}\cdot
e_{n+1}+\sum_{j=1}^{n}a_{j}y_{j}+\frac{v(y)}{w(x_{0})}$, where $a_{j}$ are
defined as before. Differentiate $u$ with respect to $y_{i}$ to obtain

$u_{i}=a_{i}+\frac{v_{i}}{w(x_{0})}\;\mbox{and}\;u_{ii}=\frac{v_{ii}}{w(x_{0})}.$

Differentiating $\log G(y)$ at $y=0$ twice yields that, for $1\leqslant
i\leqslant n$,

(3.19) $0=\frac{\beta u_{i}}{u}-\frac{(1/w)_{i}}{1/w-a}+A(y_{i}+\langle
X_{0},\epsilon_{i}\rangle)+\frac{v_{11i}}{v_{11}}$

and

(3.20) $0\geqslant\frac{\beta
u_{ii}}{u}-\beta\Big{(}\frac{u_{i}}{u}\Big{)}^{2}-\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}_{i}+A+\frac{v_{11ii}}{v_{11}}-\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-2v_{1i}^{2}-v_{ii}^{2}.$

Contracting (3.20) with $\tilde{F}^{ii}$ and using
$\tilde{F}^{ii}v_{ii}=\tilde{F}^{ii}a_{ii}=\tilde{f}$ and (3.19), similar to
(3.10), we will arrive at

(3.21) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{C\beta}{u}-\beta\tilde{F}^{ii}\Big{(}\frac{u_{i}}{u}\Big{)}^{2}-\frac{2}{v_{11}}\sum_{i\geqslant
2}\tilde{F}^{1i,i1}a_{1i1}a_{i11}-CA-Cv_{11}$
$\displaystyle\;-\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\frac{aw}{1-aw}\tilde{F}^{ii}v_{ii}^{2}+A\sum\tilde{F}^{ii}$

as long as $v_{11}$ large enough.

If $\kappa_{n}\leqslant-\delta\kappa_{1}$ for some $\delta>0$, similar to
Lemma 3.2, by the critical equation (3.19) and the Cauchy-Schwarz inequality,
we have

(3.22)
$\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}\leqslant(1+\varepsilon)\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\frac{CA^{2}}{\varepsilon}+\frac{C\beta^{2}}{\varepsilon}\frac{1}{u^{2}}$

for any $\varepsilon>0$, where $C$ depends on $|Du|_{C^{0}}$. Substituting the
above inequality into (3.21), we then have

(3.23) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{a}{2}\tilde{F}^{ii}v_{ii}^{2}+\frac{C\beta}{u}-\frac{C\beta^{2}}{\varepsilon}\frac{\sum\tilde{F}^{ii}}{u^{2}}-\frac{CA^{2}}{\varepsilon}\sum\tilde{F}^{ii}-Cv_{11}$

if $\varepsilon>0$ is sufficiently small and $v_{11}>A$. Utilizing
$\kappa_{n}\leqslant-\delta\kappa_{1}$ and Lemma 2.2 (4), we then obtain

(3.24) $\displaystyle
0\geqslant\frac{a\delta^{2}\theta}{2}v_{11}^{2}\sum\tilde{F}^{ii}+\frac{C\beta}{u}-\frac{C\beta^{2}}{\varepsilon}\frac{\sum\tilde{F}^{ii}}{u^{2}}-\frac{CA^{2}}{\varepsilon}\sum\tilde{F}^{ii}-Cv_{11},$

By Lemma 2.2 (2), we derive, for some $C$ depending on $n$, $p$, $|u|_{C^{1}}$
and $|f|_{C^{2}}$, that

$u^{2}v_{11}^{2}\leqslant C,$

which implies (3.18).

For any $\epsilon>0$ sufficiently small, choose $\delta=\delta(\epsilon)$ such
that Lemma 2.5 holds. In the following, we will assume
$\kappa_{n}\geqslant-\delta\kappa_{1}$. By Lemma 2.4 and 2.5, the inequality
(3.21) becomes

(3.25) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{C\beta}{u}-\beta\tilde{F}^{ii}\Big{(}\frac{u_{i}}{u}\Big{)}^{2}+(1+\epsilon)\sum_{i\geqslant
2}\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-CA$
$\displaystyle\;-(1+\epsilon)\sum_{i\geqslant
2}\tilde{F}^{11}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}+\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+a\tilde{F}^{ii}v_{ii}^{2}+\frac{A}{2}\sum\tilde{F}^{ii}$

as long as $A$ is sufficiently large. By (3.19) and the Cauchy-Schwarz
inequality, we have, for $i\geqslant 2$,

$-\beta\tilde{F}^{ii}\Big{(}\frac{u_{i}}{u}\Big{)}^{2}\geqslant-\frac{3}{\beta}\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-\frac{3}{\beta}\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}-\frac{CA^{2}}{\beta}\tilde{F}^{ii}.$

Taking $i=1$ and $\varepsilon=\epsilon$ in (3.22) and substituting it together
with the above inequality into (3.25), we obtain

(3.26) $\displaystyle 0\geqslant$
$\displaystyle\;a\tilde{F}^{ii}v_{ii}^{2}-\beta\tilde{F}^{11}\Big{(}\frac{u_{1}}{u}\Big{)}^{2}+(\epsilon-\frac{3}{\beta})\sum_{i\geqslant
2}\tilde{F}^{ii}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}-\frac{CA^{2}}{\beta}\sum\tilde{F}^{ii}$
$\displaystyle\;-\epsilon\tilde{F}^{11}\Big{(}\frac{(1/w)_{1}}{1/w-a}\Big{)}^{2}-\frac{CA^{2}}{\epsilon}\tilde{F}^{11}-\frac{C\beta^{2}}{\epsilon}\frac{\tilde{F}^{11}}{u^{2}}+\frac{A}{2}\sum\tilde{F}^{ii}-CA$
$\displaystyle\;-(1+\epsilon)\sum_{i\geqslant
2}\tilde{F}^{11}\Big{(}\frac{v_{11i}}{v_{11}}\Big{)}^{2}+(1-\frac{3}{\beta})\sum_{i\geqslant
2}\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}+\frac{C\beta}{u}.$

Recall that
$\tilde{F}^{11}\leqslant\tilde{F}^{22}\leqslant\cdots\leqslant\tilde{F}^{nn}$.
Choosing $\epsilon$ sufficiently small and $\beta$ sufficiently large and
using (3.22) for $i\geqslant 2$, we arrive at

(3.27) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{a}{2}\tilde{F}^{ii}v_{ii}^{2}-(3\epsilon+\frac{3}{\beta})\sum_{i\geqslant
2}\tilde{F}^{ii}\Big{(}\frac{(1/w)_{i}}{1/w-a}\Big{)}^{2}$
$\displaystyle\;+\frac{C\beta}{u}-\frac{C\beta^{2}}{\epsilon}\frac{\tilde{F}^{11}}{u^{2}}-\frac{CA^{2}}{\epsilon}\tilde{F}^{11}+\frac{A}{4}\sum\tilde{F}^{ii}-CA$
$\displaystyle\geqslant$
$\displaystyle\;\frac{a}{4}\tilde{F}^{ii}v_{ii}^{2}+\frac{C\beta}{u}-\frac{C\beta^{2}}{\epsilon}\frac{\tilde{F}^{11}}{u^{2}}-\frac{CA^{2}}{\epsilon}\tilde{F}^{11}-CA.$

By Lemma 2.2 $(1)$, we conclude

$u^{2}v_{11}^{2}\leq C$

from the last inequality of (3.27), where $C$ depends on $n$, $p$,
$|u|_{C^{1}}$, $\inf\tilde{f}$ and $|f|_{C^{2}}$. This completes the proof of
Theorem 3.3.

∎

##  4\. Boundary $C^{2}$ estimates

In this section, we establish the boundary $C^{2}$ estimate. The main theorem
in this section is stated as below.

######  Theorem 4.1.

Let $\Omega$ be a strictly convex bounded domain in $\mathbb{R}^{n}$ with
smooth boundary $\partial\Omega$. Suppose $u\in C^{3}(\bar{\Omega})$ is an
admissible solution to (1.1) and (1.2). Then, there exists a positive constant
$C$ such that

(4.1) $\max_{\partial\Omega}|D^{2}u|\leqslant C,$

where $C$ depends on $n,p$, $|u|_{C^{1}}$, $|f|_{C^{1}}$, $\inf f$ and
$\partial\Omega$.

Without loss of generality, we assume $0\in\partial\Omega$. Suppose that the
boundary $\partial\Omega$ around the origin is given by

(4.2)
$x_{n}=\rho(x^{\prime})=\frac{1}{2}\sum_{\alpha<n}\kappa_{\alpha}^{b}x_{\alpha}^{2}+O(|x^{\prime}|^{3}),$

where $x^{\prime}=(x_{1},\cdots,x_{n-1})$ and
$\kappa_{1}^{b},\cdots,\kappa_{n-1}^{b}$ are the principal curvatures of
$\partial\Omega$ at the origin. Differentiating the boundary condition $u=0$
at the origin twice, we obtain that

(4.3) $|u_{\alpha\beta}(0)|\leqslant
C,\;\mbox{for}\;1\leqslant\alpha,\beta\leqslant n-1.$

For the tangential-normal estimate, we need some idea from [26]. Let
$d(x)=dist(x,\partial\Omega)$ be the distance function to the boundary and
$\Omega_{\delta}:=\\{x\in\Omega:|x|<\delta\\}$. We consider

$W:=u_{\alpha}+\rho_{\alpha}u_{n}-\frac{1}{2}\sum_{1\leqslant\beta\leqslant
n-1}u_{\beta}^{2},$

where $1\leqslant\alpha\leqslant n-1$. Similar to Lemma 5.1 in [26] and Lemma
5.3 in [29], we can prove the following lemma.

######  Lemma 4.2.

For sufficiently small $\delta>0$, we have,

(4.4) $G^{ij}W_{ij}\leqslant
C(1+|DW|+\sum_{i}G^{ii}+G^{ij}W_{i}W_{j})\;\mbox{in}\;\Omega_{\delta},$

where $C$ depends on $n$, $p$, $|u|_{C^{1}}$, $|\tilde{f}|_{C^{1}}$,
$\inf\tilde{f}$ and $\partial\Omega$.

###### Proof.

Note that

(4.5) $\displaystyle G^{ij}W_{ij}+G^{s}W_{s}=$
$\displaystyle\;\tilde{f}_{\alpha}+\rho_{\alpha}\tilde{f}_{n}-\sum_{\beta\leqslant
n-1}u_{\beta}\tilde{f}_{\beta}+2G^{ij}u_{ni}\rho_{\alpha j}$
$\displaystyle\;-\sum_{\beta\leqslant n-1}G^{ij}u_{\beta i}u_{\beta
j}+u_{n}G^{ij}\rho_{\alpha ij}+u_{n}G^{s}\rho_{\alpha s}.$

It is easy to see that

(4.6)
$G^{ij}=\frac{1}{w}\sum_{s,t}\tilde{F}^{st}\gamma^{si}\gamma^{jt}\;\mbox{and}\;u_{ij}=w\sum_{s,t}\gamma_{is}a_{st}\gamma_{tj}.$

It follows that

$\sum_{\beta\leqslant n-1}G^{ij}u_{\beta i}u_{\beta j}=w\sum_{\beta\leqslant
n-1}\sum_{s,t}\frac{\partial\tilde{F}}{\partial a_{ij}}\gamma_{\beta
s}\gamma_{\beta t}a_{si}a_{tj}.$

Suppose $b_{ij}$ is the orthogonal matrix that diagonalize $a_{ij}$ and
$\tilde{F}^{ij}$ simultaneously, i.e.

(4.7)
$\tilde{F}^{ij}=\sum_{s}b_{is}\frac{\partial\tilde{F}}{\partial\kappa_{s}}b_{js}\;\mbox{and}\;a_{ij}=\sum_{s}b_{is}\kappa_{s}b_{js}.$

Therefore, we get

$\sum_{\beta\leqslant n-1}G^{ij}u_{\beta i}u_{\beta j}=w\sum_{\beta\leqslant
n-1}\sum_{i}\Big{(}\sum_{s}\gamma_{\beta
s}b_{si}\Big{)}^{2}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}.$

Define the matrix $\eta=(\eta_{ij})$ by

$\eta_{ij}=\sum_{s}\gamma_{is}b_{sj}.$

It is easy to verify that $\eta\cdot\eta^{T}=g$ and
$|\det(\eta)|=\sqrt{1+|Du|^{2}}$. Hence,

(4.8) $\sum_{\beta\leqslant n-1}G^{ij}u_{\beta i}u_{\beta
j}=w\sum_{\beta\leqslant n-1}\sum_{i}(\eta_{\beta
i})^{2}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}.$

By (4.6) and (4.7), we obtain

(4.9) $|G^{ij}u_{ni}\rho_{\alpha j}|\leqslant
C\sum_{i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}|\kappa_{i}|\;\mbox{and}\;\tilde{F}^{ij}a_{it}=\sum_{i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}b_{ji}b_{ti}.$

By (2.5) and the above equality, we derive that

$G^{s}\rho_{\alpha s}\leqslant
C\sum_{i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}|\kappa_{i}|.$

A direct calculation shows that

$|\tilde{f}_{\alpha}+\tilde{f}_{n}\rho_{\alpha}-\sum_{\beta\leqslant
n-1}u_{\beta}\tilde{f}_{\beta}|\leqslant C(1+|DW|).$

We therefore arrive at

(4.10) $\displaystyle G^{ij}W_{ij}\leqslant$
$\displaystyle\;C\Big{(}1+|DW|+\sum G^{ii}\Big{)}-G^{s}W_{s}$
$\displaystyle\;+C\sum_{i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}|\kappa_{i}|-\sum_{\beta\leqslant
n-1}G^{ij}u_{\beta i}u_{\beta j}.$

Next we estimate $G^{s}W_{s}$. By Lemma 2.1 and (4.9), we see that

(4.11) $\displaystyle-G^{s}W_{s}=$
$\displaystyle\;\frac{u_{s}W_{s}}{w^{2}}\tilde{f}+\frac{2}{w}\sum_{t,i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}(b_{ti}u_{t})(b_{ji}\gamma^{sj})W_{s}$
$\displaystyle\leqslant$
$\displaystyle\;C|DW|+\frac{2}{w}\sum_{t,i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}(b_{ti}u_{t})(b_{ji}\gamma^{sj})W_{s}.$

Now we divide the proof into two cases: (a) $\sum_{\beta\leqslant
n-1}\eta_{\beta i}^{2}\geqslant\epsilon$ for all index $i$; (b)
$\sum_{\beta\leqslant n-1}\eta_{\beta r}^{2}<\epsilon$ for some index
$1\leqslant r\leqslant n$, where $\epsilon>0$ is a small constant to be
determined later.

Case (a). By our assumption and (4.8), we have

(4.12) $\sum_{\beta\leqslant n-1}G^{ij}u_{\beta i}u_{\beta
j}\geqslant\epsilon\sum_{i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}.$

The second term in the right hand side of the inequality (4.11) can be
estimated as

(4.13)
$\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}(b_{ti}u_{t})(b_{ji}\gamma^{sj})W_{s}\leqslant\frac{\epsilon}{4}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}+\frac{C}{\epsilon}\frac{\partial\tilde{F}}{\partial\kappa_{i}}(b_{ji}\gamma^{sj}W_{s})^{2}$

by the Cauchy-Schwarz inequality. Note that

(4.14)
$\frac{\partial\tilde{F}}{\partial\kappa_{i}}(b_{ji}\gamma^{sj}W_{s})^{2}=G^{ij}W_{i}W_{j}.$

Combining (4.14), (4.13) with (4.11), we obtain

(4.15) $-G^{s}W_{s}\leqslant
C(|DW|)+\frac{\epsilon}{2}\sum_{i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}+\frac{C}{\epsilon}G^{ij}W_{i}W_{j}.$

Substituting (4.15) and (4.12) into (4.10) and applying the Cauchy-Schwarz
inequality to $\sum\frac{\partial\tilde{F}}{\partial\kappa_{i}}|\kappa_{i}|$,
we get the desired inequality (4.4).

Case (b). We may assume that the principal curvatures are ordered as
$\kappa_{1}\geqslant\kappa_{2}\geqslant\cdots\geqslant\kappa_{n}$. Then, we
have
$\frac{\partial\tilde{F}}{\partial\kappa_{1}}\leqslant\frac{\partial\tilde{F}}{\partial\kappa_{2}}\leqslant\cdots\leqslant\frac{\partial\tilde{F}}{\partial\kappa_{n}}$.
Taking into account the fact that
$|\eta_{nr}|\leqslant\sqrt{1+|Du|_{C^{0}}^{2}}$, we see that

$1\leqslant|\det(\eta)|\leqslant\sqrt{1+|Du|_{C^{0}}^{2}}|\det(\eta^{\prime})|+C\epsilon,$

where $\eta^{\prime}=\\{\eta_{pq}\\}_{p\neq n,q\neq r}$, which implies

$|\det(\eta^{\prime})|\geqslant\frac{1}{2\sqrt{1+|Du|_{C^{0}}^{2}}}$

for sufficiently small $\epsilon$. On the other hand, for any fixed index
$i\neq r$, it holds that

$|\det(\eta^{\prime})|\leqslant C\sum_{\beta\neq n}|\eta_{\beta i}|$

for some positive constant $C$ depending on $n$ and $|Du|_{C^{0}}$. By the
above two inequalities, we derive that, for any $i\neq r$,

$\sum_{\beta\leqslant n-1}\eta_{\beta i}^{2}\geqslant c_{1}$

for some positive constant $c_{1}$ depending on $n$ and $|Du|_{C^{0}}$.
Substituting the above inequality into (4.8), we obtain

(4.16) $\sum_{\beta\leqslant n-1}G^{ij}u_{\beta i}u_{\beta j}\geqslant
c_{1}\sum_{i\neq
r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}.$

If $\kappa_{r}\leqslant 0$, by the same argument as that of Lemma 2.20 in [15]
and $\sum\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}\geqslant 0$,
we obtain

$\sum_{i\neq
r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}\geqslant\frac{1}{n}\sum_{i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}.$

We then get the desired inequality (4.4) by the same argument as Case (a) with
$\epsilon$ in (4.12) replaced by $\frac{c_{1}}{n}$. Therefore, we may assume
$\kappa_{r}>0$. Note that $\kappa_{r}\geqslant\kappa_{n}$. It is easy to see
that

$(b_{jr}\gamma^{sj})W_{s}=w\Big{(}\eta_{\alpha
r}+\rho_{\alpha}\eta_{nr}-\sum_{\beta\leq n-1}u_{\beta}\eta_{\beta
r}\Big{)}\kappa_{r}+b_{jr}\gamma^{sj}\rho_{\alpha s}u_{n},$

which implies that

(4.17) $(b_{jr}\gamma^{sj})W_{s}\leq
Cw(\epsilon+|\rho_{\alpha}|)\kappa_{r}+C.$

First, suppose that $\kappa_{n}\geqslant-\frac{\kappa_{r}}{2(p-1)}$. We see
$\kappa_{r}+\kappa_{i_{1}}+\cdots+\kappa_{i_{p-1}}\geqslant\frac{\kappa_{r}}{2}$.
It follows that

(4.18) $\displaystyle\frac{\partial F}{\partial\kappa_{r}}=$
$\displaystyle\;\sum_{r\notin\\{i_{1},\cdots,i_{p-1}\\}}\frac{F(\kappa)}{\kappa_{r}+\kappa_{i_{1}}+\cdots+\kappa_{p-1}}\leqslant\frac{C}{\kappa_{r}}.$

By (4.18) and the Cauchy-Schwarz inequality, we have

$\displaystyle\;\sum_{t,i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}(b_{ti}u_{t})(b_{ji}\gamma^{sj})W_{s}\leqslant
C|DW|+\frac{\epsilon}{4}\sum_{i\neq
r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}+\frac{C}{\epsilon}G^{ij}W_{i}W_{j}.$

Substituting the above inequality into (4.11), we obtain

(4.19) $-G^{s}W_{s}\leqslant C|DW|+\frac{\epsilon}{4}\sum_{i\neq
r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}+\frac{C}{\epsilon}G^{ij}W_{i}W_{j}.$

By (4.19), (4.18), (4.16) and applying the Cauchy-Schwarz inequality to
$\sum_{i\neq r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}|\kappa_{i}|$, we
conclude (4.4) from (4.10).

Now suppose that $\kappa_{n}\leqslant-\frac{\kappa_{r}}{2(p-1)}$. Note that
$r\neq n$ in this case. Since
$\frac{\partial\tilde{F}}{\partial\kappa_{r}}\kappa_{r}=\tilde{f}-\sum_{i\neq
r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}$, we have, by (4.17)
and the Cauchy-Schwarz inequality,

$\displaystyle\;\sum_{t,i}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}(b_{ti}u_{t})(b_{ji}\gamma^{sj})W_{s}$
$\displaystyle\leqslant$
$\displaystyle\;C\frac{\partial\tilde{F}}{\partial\kappa_{r}}\kappa_{r}+C(\epsilon+|\rho_{\alpha}|)\frac{\partial\tilde{F}}{\partial\kappa_{r}}\kappa_{r}^{2}+\frac{\epsilon}{4}\sum_{i\neq
r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}+\frac{C}{\epsilon}G^{ij}W_{i}W_{j}$
$\displaystyle\leqslant$ $\displaystyle\;C+\frac{C}{\epsilon}\sum_{i\neq
r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}+C(\epsilon+|\rho_{\alpha}|)\frac{\partial\tilde{F}}{\partial\kappa_{n}}\kappa_{n}^{2}+\frac{\epsilon}{2}\sum_{i\neq
r}\frac{\partial\tilde{F}}{\partial\kappa_{i}}\kappa_{i}^{2}+\frac{C}{\epsilon}G^{ij}W_{i}W_{j}.$

Combining (4.10), (4.11) and (4.16) with the above inequality and choosing
$\epsilon,\delta$ sufficiently small such that $C(\epsilon+\delta)<c_{1}$, we
can get (4.4). Thus, the proof of Lemma 4.2 is finished.

∎

With the above lemma in hand, we can construct a suitable barrier function to
prove the following estimate

(4.20) $|u_{\alpha n}(0)|\leqslant C,\;\mbox{for}\;1\leqslant\alpha\leqslant
n-1.$

Since the boundary $\partial\Omega$ is smooth and strictly convex, we can find
a smooth strictly convex function $v$ satisfying $v=0$ on $\partial\Omega$ and
$v\leqslant 0$ in $\Omega$. Therefore, we have $D^{2}v(x)\in P_{1}$ for
$x\in\overline{\Omega}$. Moreover, we can find a small positive constant
$\theta$ such that $\tilde{v}=v-\theta|x|^{2}$ satisfies $D^{2}\tilde{v}(x)\in
P_{1}$ for $x\in\overline{\Omega}$. Now we can finish the proof of (4.20).
Consider the following barrier function

(4.21) $\Psi=\tilde{v}+\frac{\theta}{2}|x|^{2}-td+\frac{N}{2}d^{2},$

where $t$ and $N$ are two positive constants to be chosen later. A direct
calculation shows that

(4.22) $\displaystyle G^{ij}\Psi_{ij}=G^{ij}(D^{2}\tilde{v}+NDd\otimes
Dd)_{ij}+G^{ij}(Nd-t)d_{ij}+\theta\sum G^{ii}.$

Concerning $D^{2}\tilde{v}\in P_{1}$, we see that $D^{2}\tilde{v}+NDd\otimes
Dd\in P_{1}$. By the concavity and homogeneity of $G$ with respect to
$D^{2}u$, we see that

(4.23) $G^{ij}(D^{2}\tilde{v}+NDd\otimes Dd)_{ij}\geq
G(D^{2}\tilde{v}+NDd\otimes Dd,Du)=\tilde{F}(V_{u}),$

where $V_{u}:=\frac{1}{w}\\{\gamma^{ik}V_{kl}\gamma^{lj}\\}$ and
$V:=\\{\tilde{v}_{kl}+Nd_{k}d_{l}\\}$. We adopt the convention that the
eigenvalues are arranged in algebraically nondecreasing order:

$\lambda_{1}\leqslant\lambda_{2}\leqslant\cdots\leqslant\lambda_{n}.$

Then, by the Weyl theorem, it is easy to see that
$\lambda_{\alpha}(V)\geqslant\lambda_{\alpha}(D^{2}\tilde{v})$ for
$\alpha=1,2,\cdots,n-1$ and
$\lambda_{n}(V)\geqslant\lambda_{1}(D^{2}\tilde{v})+N$. Denote
$\theta_{0}=\inf_{x\in\Omega}\frac{\lambda_{1}(D^{2}\tilde{v})}{\lambda_{n}(D^{2}\tilde{v})}$.
We further derive that
$\lambda_{\alpha}(V)\geqslant\theta_{0}\lambda_{\alpha}(D^{2}\tilde{v})$ for
$\alpha=1,2,\cdots,n-1$ and
$\lambda_{n}(V)\geqslant\theta_{0}(\lambda_{n}(D^{2}\tilde{v})+N)$. On the
other hand, by the Ostrowski Theorem, there exist $\theta_{k}>0$ such that
$\lambda_{k}(V_{u})=\theta_{k}\lambda_{k}(V)$ for $k=1,2,\cdots,n$, where
$\theta_{k}$ is uniformly positive only depending on $|Du|_{C^{0}}$.
Substituting (4.23) into (4.22) and by Lemma 2.3, we have, for sufficiently
large $N$ and sufficiently small $\delta$ and $t$,

$G^{ij}\Psi_{ij}\geqslant\bar{\theta}N^{\frac{p}{n}}+\frac{\theta}{2}\sum
G^{ii}\geqslant N^{\frac{p-1}{n}}+\frac{\theta}{2}\sum G^{ii},$

where $\bar{\theta}$ depends on $\theta_{0},\theta_{1},\cdots\theta_{n}$ and
$\lambda(D^{2}\tilde{v})$. Also, note that $|D\Psi|=|D\tilde{v}+\theta
x-tDd+NdDd|\leqslant C+N\delta$ in $\Omega_{\delta}$. Choosing $\delta$
sufficiently small, we obtain

$G^{ij}\Psi_{ij}\geqslant N^{\frac{1}{n}}|D\Psi|+\frac{\theta}{2}\Big{(}\sum
G^{ii}+1\Big{)}.$

Let $\tilde{W}=1-e^{-BW}$, where constant $B>0$ is sufficiently large. It
follows that

(4.24)
$G^{ij}(R\Psi-\tilde{W})_{ij}\geqslant\sqrt[n]{N}(|DR\Psi|-|D\tilde{W}|)$

if we choose $R\gg B\gg 1$. Since
$\Psi|_{\partial\Omega_{\delta}\cap\partial\Omega}\leqslant-\frac{\theta}{2}|x^{\prime}|^{2}$,
$\Psi|_{\partial\Omega_{\delta}\cap\Omega}\leqslant-\frac{\theta}{2}\delta^{2}$
and
$-W|_{\partial\Omega}=\frac{1}{2}\sum_{\beta}u_{\beta}^{2}=\frac{1}{2}\sum_{\beta}\rho_{\beta}^{2}u_{n}^{2}\leqslant
C|x^{\prime}|^{2}$ for some $C$ depending on $|\rho|_{C^{2}}$ and
$|u|_{C^{1}}$, it is easy to verify that on $\partial\Omega_{\delta}$ we have
$R\Psi-\tilde{W}\leqslant 0$ if $R$ is large enough. By the maximum principle
we get

$R\Psi-\tilde{W}\leqslant 0\;\mbox{in}\;\Omega_{\delta}.$

As $(R\Psi-\tilde{W})(0)=0$ we deduce that $u_{\alpha n}(0)\geqslant-C$.
Consider
$W^{\prime}=-u_{\alpha}-\rho_{\alpha}u_{n}-\frac{1}{2}\sum_{1\leqslant\beta\leqslant
n-1}u_{\beta}^{2}$. By the same argument we deduce $u_{\alpha n}(0)\leqslant
C$. The proof of (4.20) is completed.

Finally, we prove the double normal estimate

(4.25) $|u_{nn}(0)|\leqslant C.$

We only need to derive an upper bound for $u_{nn}(0)$ by
$\sum_{i}\kappa_{i}>0$ since $P_{p}\subset P_{n}$. Actually, by the inequality
of arithmetic and geometric means, we have
$\frac{p\sum_{i}\kappa_{i}}{n}\geqslant\tilde{F}(\kappa)\geqslant\inf\tilde{f}>0$.
By (1.8) in [6], i.e. $ad(x)\leqslant-u(x)\leqslant\frac{1}{a}d(x)$ in
$\Omega$ for some constant $a>0$ depending on $\inf\tilde{f}$, we obtain

$u_{n}(0)\leqslant-a.$

We have, with respect to a principal coordinate system at the origin,
$u_{\alpha\beta}(0)=-u_{n}\kappa_{\alpha}^{b}\delta_{\alpha\beta}$ for
$1\leqslant\alpha,\beta\leqslant n-1$ and
$g^{ij}=\delta_{ij}-\frac{|Du|^{2}}{w^{2}}\delta_{in}\delta_{jn}$. Therefore,
the matrix in (2.2) has the following form

$\left(\begin{array}[]{ccccc}-u_{n}\kappa_{1}^{b}&0&\cdots&0&\frac{u_{1n}}{w}\\\
0&-u_{n}\kappa_{2}^{b}&\cdots&0&\frac{u_{2n}}{w}\\\
\vdots&\vdots&\ddots&\vdots&\vdots\\\
0&0&\cdots&-u_{n}\kappa_{n-1}^{b}&\frac{u_{n-1n}}{w}\\\
\frac{u_{n1}}{w}&\frac{u_{n2}}{w}&\cdots&\frac{u_{nn-1}}{w}&\frac{u_{nn}}{w^{2}}\end{array}\right).$

By Lemma 1.2 of [4], when $|u_{nn}|$ goes to infinity, the eigenvalues
$w(\kappa_{1},\cdots,\kappa_{n})$ of the above matrix behave like

$w\kappa_{\alpha}=-u_{n}\kappa_{\alpha}^{b}+o(1),\;1\leqslant\alpha\leqslant
n-1,\;w\kappa_{n}=\frac{u_{nn}}{w^{2}}+O(1),$

where $o(1)$ and $O(1)$ are uniform only depending on
$-u_{n}\kappa^{b}_{\alpha}$ and $|u_{n\alpha}(0)|$. Since $-u_{n}(0)\geqslant
a$ and $(\kappa_{1}^{b},\cdots,\kappa_{n-1}^{b})\in\mathcal{P}_{p}$, there
exists a uniform positive constant $\varsigma$ such that

$\Pi_{1\leqslant i_{1}<\cdots<i_{p}\leqslant
n-1}(\kappa_{i_{1}}+\cdots+\kappa_{i_{p}})\geqslant\varsigma,$

as long as $u_{nn}$ is large enough. Returning to the equation (1.1), we have

$\Pi_{1\leqslant i_{1}<\cdots<i_{p-1}\leqslant
n-1}(\kappa_{n}+\kappa_{i_{1}}+\cdots+\kappa_{i_{p-1}})\leqslant\frac{C}{\varsigma},$

from which we can derive an upper bound for $\kappa_{n}$. Therefore,
$u_{nn}(0)$ is bounded from above. Combining (4.3), (4.20) with (4.25), we
obtain (4.1).

##  5\. Gradient Estimates

In this section, we establish the gradient estimate. First, by the existence
of a subsolution $\underline{u}$ and the maximum principle, it is easy to
derive

(5.1) $\sup_{\Omega}|u|+\sup_{\partial\Omega}|Du|\leqslant C,$

where $C$ depends on $|\underline{u}|_{C^{1}}$. Next, we prove the global
gradient estimate following the argument in [6]. We have

######  Theorem 5.1.

Let $u\in C^{3}(\Omega)\cap C^{1}(\bar{\Omega})$ be an admissible solution to
the problem (1.1) and (1.2). Suppose $f(x,z,\nu)\in
C^{\infty}(\bar{\Omega}\times\mathbb{R}\times\mathbb{S}^{n})>0$ and
$f_{z}\geqslant 0.$ Then, the estimate

(5.2) $\sup_{\Omega}|Du|\leqslant C(1+\sup_{\partial\Omega}|Du|)$

holds for a positive constant $C$ depending on $n$, $p$, $|u|_{C^{0}}$, $\inf
f$ and $|f|_{C^{1}}$.

###### Proof.

As in [6], we consider the following quantity

$Q=|Du|e^{Au},$

where $A$ is a large constant to be determined later. Suppose that $Q$
achieves its maximum at $x_{0}\in\Omega$. At this point we may rotate the
coordinates such that

$u_{1}=|Du|>0\;\mbox{ and}\;u_{\alpha}=0\;\mbox{for}\;\alpha\geqslant 2.$

Differentiating $\log Q$ at $x_{0}$ once, we obtain

(5.3) $\frac{u_{1i}}{u_{1}}+Au_{i}=0\;\mbox{for}\;i=1,\cdots,n.$

By the above equality, we may further rotate the coordinates so that $u_{ij}$
is diagonal at $x_{0}$. From (2.2), we find that

$a_{11}=\frac{u_{11}}{w^{3}},\;a_{ii}=\frac{u_{ii}}{w}\;\mbox{for}\;i\geqslant
2\;\mbox{and}\;a_{ij}=0\;\mbox{for}\;i\neq j.$

Denote by $\alpha_{ij}$ the positive semidefinite matrix defined as

$\alpha_{11}=\frac{1}{w^{3}},\;\alpha_{1i}=\frac{1}{w^{2}},\;i\geqslant
2,\;\alpha_{ij}=\frac{1}{w},\;j\geqslant i\geqslant 2.$

As in [49], denote by $\hat{F}^{ij}=F^{ij}\alpha_{ij}$ the Hadamard product of
the matrix $F^{ij}$ and $\alpha_{ij}$, which is positive definite since
$F^{ij}$ is diagonal and positive definite. Differentiating $\log Q$ at
$x_{0}$ a second time and contracting with $\hat{F}^{ij}$, we get

(5.4)
$0\geqslant\frac{\hat{F}^{ij}u_{1ij}}{u_{1}}-\frac{\hat{F}^{ij}u_{1i}u_{1j}}{u_{1}^{2}}+A\hat{F}^{ij}u_{ij}.$

Differentiating $a_{ij}$ and the equation (1.1), we obtain

(5.5)
$F^{ij}a_{ij1}=\hat{F}^{ij}u_{ij1}-\frac{2F^{11}u_{11}^{2}u_{1}}{w^{5}}-\frac{u_{11}u_{1}}{w^{2}}F^{ij}a_{ij}=(f)_{1}.$

Note that $(f)_{1}=f_{x_{1}}+f_{u}u_{1}+f_{\nu_{j}}(\nu_{j})_{1}$. By (5.3)
and $f_{u}\geqslant 0$, we have

(5.6)
$\frac{(f)_{1}}{u_{1}}\geqslant-\frac{C}{u_{1}}+\frac{u_{11}f_{\nu_{1}}}{u_{1}w}\Big{(}\frac{u_{1}^{2}}{w^{2}}-1\Big{)}-f_{\nu_{n+1}}\frac{u_{11}}{w^{3}}\geqslant-\frac{CA}{u_{1}},$

where $C$ depends on
$|Df|_{C^{0}}\equiv\sup_{\bar{\Omega}\times[-C_{0},C_{0}]\times\mathbb{S}^{n}}|Df|$.
Substituting (5.5) and (5.6) into (5.4), we have

(5.7) $\displaystyle 0\geqslant$
$\displaystyle\;\frac{2F^{11}u_{11}^{2}}{w^{5}}+\frac{u_{11}}{w^{2}}F^{ij}a_{ij}-\frac{CA}{u_{1}}-\frac{\hat{F}^{11}u_{11}^{2}}{u_{1}^{2}}+A\hat{F}^{ij}u_{ij}$
$\displaystyle\geqslant$
$\displaystyle\;(\frac{2u_{1}^{2}}{w^{2}}-1)\frac{F^{11}u_{11}^{2}}{w^{3}u_{1}^{2}}+(1-\frac{u_{1}^{2}}{w^{2}})AC_{n}^{p}f-\frac{CA}{u_{1}},$

where in the last inequality we used
$F^{ij}a_{ij}=\hat{F}^{ij}u_{ij}=C_{n}^{p}f$. From (5.3), we see
$u_{11}=-Au_{1}^{2}<0$. By Lemma 2.2 (4) and (2), there exists a positive
constant $c_{2}$ depending on $n$, $p$ and $\inf f$ such that

(5.8) $0\geqslant\frac{c_{2}A^{2}u_{1}^{4}}{2w^{5}}-\frac{CA}{u_{1}},$

which is a contradiction to large $A$ and $|Du|$. Therefore, $Q$ attains its
maximum on $\partial\Omega$ and hence (5.2) is proved.

∎

######  Remark 5.2.

The interior gradient estimate was proved by Li [32] for Weingarten curvature
equations of general type by the normal perturbation method.

## References

  * [1] I. Bakelman and B. Kantor, Existence of a hypersurface homeomorphic to the sphere in Euclidean space with a given mean curvature. (Russian) Geometry and topology, No. 1 (Russian), pp. 3-10. Leningrad. Gos. Ped. Inst. im. Gercena, Leningrad, 1974.
  * [2] J. Ball, Differentiability properties of symmetric and isotropic functions, Duke Math. J. 51 (1984), 699-728.
  * [3] L. Caffarelli, L. Nirenberg and J. Spruck, The Dirichlet problem for nonlinear second order elliptic equations I. Monge-Ampère equations, Comm. Pure Appl. Math. 37 (1984), 369-402.
  * [4] L. Caffarelli, L. Nirenberg and J. Spruck, Dirichlet problem for nonlinear second order elliptic equations III, Functions of the eigenvalues of the Hessian, Acta Math. 155 (1985), 261-301.
  * [5] L. Caffarelli, L. Nirenberg and J. Spruck, Nonlinear second order elliptic equations IV: Starshaped compact Weigarten hypersurfaces, Current topics in partial differential equations, Y.Ohya, K.Kasahara and N.Shimakura (eds), Kinokunize, Tokyo, 1985, 1-26.
  * [6] L. Caffarelli, L. Nirenberg and J. Spruck, Nonlinear second-order elliptic equations. V. The Dirichlet problem for Weingarten hypersurfaces, Comm. Pure Appl. Math. 41 (1988), no. 1, 47–70.
  * [7] D. Chen, H. Li and Z. Wang, Starshaped compact hypersurfaces with prescribed Weingarten curvature in warped product manifolds, Calc. Var. Partial Differential Equations 57 (2018), Paper No. 42, 26 pp.
  * [8] J. Chu and H. Jiao, Curvature estimates for a class of Hessian type equations, Calc. Var. Partial Differential Equations 60 (2021), Paper No. 90, 18 pp.
  * [9] S. Dinew, Interior estimates for $p$-plurisubharmonic functions, arXiv:2006.12979.
  * [10] W. Dong, Curvature estimates for $p$-convex hypersurfaces of prescribed curvature, arXiv:2107,13387v2, to appear in Rev. Mat. Iberoam. https://doi.org/10.4171/RMI/1348.
  * [11] J. Fu, Z. Wang and D. Wu, Form-type Calabi-Yau equations, Math. Res. Lett. 17 (2010), 887-903.
  * [12] J. Fu, Z. Wang and D. Wu, Form-type equations on Kähler manifolds of nonnegative orthogonal bisectional curvature, Calc. Var. Partial Differential Equations 52 (2015), 327-344.
  * [13] C. Gerhardt, Hypersurfaces of prescribed Weingarten curvature, Math. Z. 224 (1997), 167-194.
  * [14] D. Gilbarg and N. Trudinger, Elliptic partial differential equations of second order. Second edition. Grundlehren der Mathematischen Wissenschaften [Fundamental Principles of Mathematical Sciences], 224. Springer-Verlag, Berlin, 1983. xiii+513 pp. ISBN: 3-540-13025-X.
  * [15] B. Guan, Second-order estimates and regularity for fully nonlinear elliptic equations on Riemannian manifolds, Duke Math. J. 163 (2014), 1491-1524.
  * [16] B. Guan and P. Guan, Convex Hypersurfaces of Prescribed Curvature, Ann. of Math. (2) 156 (2002), 655-674.
  * [17] B. Guan and H. Jiao, Second order estimates for Hessian type fully nonlinear elliptic equations on Riemannian manifolds, Calc. Var. Partial Differential Equations 54 (2015), 2693-2712.
  * [18] B. Guan and J. Spruck, Locally convex hypersurfaces of constant curvature with boundary. Comm. Pure Appl. Math. 57 (2004), no. 10, 1311-1331.
  * [19] P. Guan, J. Li and Y. Li, Hypersurfaces of prescribed curvature measure, Duke Math. J. 161 (2012), 1927-1942.
  * [20] P. Guan, C. Lin and X.-N. Ma, The existence of convex body with prescribed curvature measures, Int. Math. Res. Not. IMRN 2009, 1947-1975.
  * [21] P. Guan, C. Ren and Z. Wang, Global $C^{2}$-estimates for convex solutions of curvature equations, Comm. Pure Appl. Math. 68 (2015), 1287-1325.
  * [22] P. Guan and X. Zhang, A class of curvature equations, Pure Appl. Math. Q. 17 (2021), no. 3, 865-907.
  * [23] F. Harvey and H. Lawson, An introduction to potential theory in calibrated geometry, Amer. J. Math. 131 (2009), no. 4, 893-944.
  * [24] F. Harvey and H. Lawson, Geometric plurisubharmonicity and convexity: an introduction, Adv. Math. 230 (2012), 2428-2456.
  * [25] F. Harvey and H. Lawson, p-convexity, p-plurisubharmonicity and the Levi problem, Indiana Univ. Math. J. 62 (2013), 149-169.
  * [26] N. Ivochkina, Solution of the Dirichlet problem for equations of $m$th order curvature. (Russian), Mat. Sb. 180 (1989), 867-887, 991; translation in Math. USSR-Sb. 67 (1990), no. 2, 317-339.
  * [27] N. Ivochkina, The Dirichlet problem for the curvature equation of order $m$, Algebra i Analiz 2 (1990), 192-217; translation in Leningrad Math. J. 2 (1991), no. 3, 631-654.
  * [28] N. Ivochkina, M. Lin and N. Trudinger, The Dirichlet problem for the prescribed curvature quotient equations with general boundary values. Geometric analysis and the calculus of variations, 125–141, Int. Press, Cambridge, MA, 1996.
  * [29] H. Jiao and Z. Sun, The Dirichlet problem for a class of prescribed curvature equations, arXiv:2202.05003v1.
  * [30] Q. Jin and Y. Li, Starshaped compact hypersurfaces with prescribed $k$-th mean curvature in hyperbolic space, Discrete Contin. Dyn. Syst. 15 (2006), no. 2, 367-377.
  * [31] A. Li, Spacelike hypersurfaces with constant Gauss-Kronecker curvature in the Minkowski space, Arch. Math. (Basel) 64 (1995), no. 6, 534-551.
  * [32] Y. Li, Interior gradient estimates for solutions of certain fully nonlinear elliptic equations, J. Differential Equations 90 (1991), 172-185.
  * [33] Y. Li and V. Oliker, Star shaped compact hypersurfaces with prescribed $m$-th mean curvature in elliptic space, J. Partial Differ. Equ. 15 (2002), 68-80.
  * [34] M. Lin and N. Trudinger, The Dirichlet problem for the prescribed curvature quotient equations, Topol. Methods Nonlinear Anal. 3 (1994), no. 2, 307-323.
  * [35] P.-L., Lions, Sur les équations de Monge-Ampère. (French) [On Monge-Ampère equations] Arch. Rational Mech. Anal. 89 (1985), no. 2, 93-122.
  * [36] V. Oliker, Hypersurfaces in $\mathbf{R}^{n+1}$ with prescribed Gaussian curvature and related equations of Monge-Ampère type, Commun. Partial Diff. Equations 9 (1984), 807-838.
  * [37] C. Ren and Z. Wang, On the curvature estimates for Hessian equations, Amer. J. Math. 141 (2019), 1281-1315.
  * [38] C. Ren and Z. Wang, The global curvature estimate for the $n-2$ Hessian equation, arXiv:2002.08702.
  * [39] C. Ren and Z. Wang, Notes on the curvature estimates for Hessian equations, arXiv:2003.14234v1.
  * [40] C. Ren, Z. Wang and L. Xiao, The prescribed curvature problem for entire hypersurfaces in Minkowski space, arXiv:2007.04493v1.
  * [41] J. Sha, $p$-convex Riemannian manifolds, Invent. Math. 83 (1986), 437-447.
  * [42] J. Sha, Handlebodies and p-convexity, J. Differential Geom. 25 (1987), 353-361.
  * [43] W. Sheng, J. Urbas and X.-J. Wang, Interior curvature bounds for a class of curvature equations, Duke Math. J. 123 (2004), no. 2, 235-264.
  * [44] G. Székelyhidi, V. Tosatti and B. Weinkove, Gauduchon metrics with prescribed volume form, Acta Math. 219 (2017), 181-211.
  * [45] V. Tosatti and B. Weinkove, The Monge-Ampère equation for $(n-1)$-plurisubharmonic functions on a compact Kähler manifold, J. Amer. Math. Soc. 30 (2017), 311-346.
  * [46] V. Tosatti and B. Weinkove, Hermitian metrics, $(n-1,n-1)$ forms and Monge-Ampère equations, J. Reine Angew. Math. 755 (2019), 67-101.
  * [47] A. Treibergs and S. Wei, Embedded hypersurfaces with prescribed mean curvature, J. Differential Geom. 18 (1983), 513-521.
  * [48] N. Trudinger, The Dirichlet problem for the prescribed curvature equations, Arch. Rational Mech. Anal., 111 (1990), 153-179.
  * [49] X.-J. Wang, Interior gradient estimates for mean curvature equations, Math. Z. 228 (1998), no. 1, 73-81.
  * [50] Z. Wang and L. Xiao, Entire spacelike hypersurfaces with constant $\sigma_{k}$ curvature in Minkowski space, arXiv:2007.01495v1.
  * [51] H. Wu, Manifolds of partially positive curvature, Indiana Univ. Math. J. 36 (1987), 525-548.


# Heterogeneity Learning for SIRS model: an Application to the COVID-19

Guanyu Hu Junxian Geng

###### Abstract

We propose a Bayesian Heterogeneity Learning approach for Susceptible-
Infected-Removal-Susceptible (SIRS) model that allows underlying clustering
patterns for transmission rate, recovery rate, and loss of immunity rate for
the latest coronavirus (COVID-19) among different regions. Our proposed method
provides simultaneously inference on parameter estimation and clustering
information which contains both number of clusters and cluster configurations.
Specifically, our key idea is to formulates the SIRS model into a hierarchical
form and assign the Mixture of Finite mixtures priors for heterogeneity
learning. The properties of the proposed models are examined and a Markov
chain Monte Carlo sampling algorithm is used to sample from the posterior
distribution. Extensive simulation studies are carried out to examine
empirical performance of the proposed methods. We further apply the proposed
methodology to analyze the state level COVID-19 data in U.S.

keywords: Bayesian Nonparametric, Cluster Learning, Infectious Diseases, MCMC,
Mixture of Finite Mixtures

##  1 Introduction

The Coronavirus Disease 2019 (COVID-19) has created a profound public health
emergency around world. It has become an epidemic with more than 5,000,000
confirmed infections worldwide as on May 21 2020. The spreading speed of
COVID-19 which is caused by a new coronavirus is faster than severe acute
respiratory syndrome (SARS) and Middle East respiratory syndrome (MERS).
Recently, the risk of COVID-19 has been a significant public-health concern
and people pay more attention on precise and timely estimates and predictions
of COVID-19. The Susceptible-Infectious-Recovered (SIR) model and its
variation approaches, such as Susceptible-Infected-Removal-Susceptible (SIRS)
(Kermack and McKendrick, 1932, 1933) and Susceptible-Exposed-Infected-Removal
(SEIR) model (Hethcote, 2000), have been widely discussed to study the
dynamical evolution of an infectious disease in a certain region. There are
rich literatures producing early results on COVID-19 based on SIR model and
its variations (Wu et al., 2020; Read et al., 2020; Tang et al., 2020). From
statistician’s perspectives, building a time-varying model under SIR and its
variations is also fully discussed for COVID-19 (Chen et al., 2020; Sun et
al., 2020; Jo et al., 2020). In most existing literature, people focus more on
dynamic regimes of the SIR models for COVID-19. They lack discussions on
heterogeneity pattern of COVID-19 among different regions.

The aim of this paper is to propose a new hierarchical SIRS model for
detecting heterogeneity pattern of COVID-19 among different regions under a
Bayesian framework. Bayesian nonparametric methods such as Dirichlet process
(DP) offer choices to do simultaneously inference on parameters’ estimation
and parameters’ heterogeneity information which contains the number of
clusters and clustering configurations. Compared with existing approaches such
as finite mixtures models, Bayesian nonparametric approach does not need to
pre-specify the number of clusters, which provides probabilistic framework for
simultaneous inference of the number of clusters and the clustering labels.
Miller and Harrison (2013) points out that the estimation of the number of
clusters under Dirichlet process mixture (DPM) model is inconsistent which
will produce extremely small clusters. One remedy for over-clustering problem
under DPM is mixture of finite mixtures model (MFM) (Miller and Harrison,
2018). The clustering properties of MFM are fully discussed in Miller and
Harrison (2018); Geng et al. (2019) and it has been widely applied in
different areas such as regional economics (Hu et al., 2020), environmental
science (Geng et al., 2019), and social science (Geng et al., 2019). Thus, the
key idea of this paper is to assign MFM priors on different parameters of the
SIRS model to capture the heterogeneity of parameters among different regions.
The contribution of this paper are two-fold. First, we formulate a Bayesian
heterogeneity learning model for SIRS under MFM. To our best knowledge, this
is the first time when MFM is applied into epidemiology models such as SIRS.
Our proposed Bayesian approach successfully captures the heterogeneity of
three different parameters under the SIRS model among different regions while
also considering uncertainty in estimation of the number of clusters. Several
interesting findings based on our proposed method are discovered for COVID-19
data in US.

This paper is organized as follows. Section 2 presents the motivating data we
analyze. We discuss our proposed Bayesian hierarchical model for heterogeneity
learning under SIRS model framework in Section 3. The performance of our
proposed method is illustrated via simulation studies in Section 4. Section 5
is devoted to the analysis of state level COVID-19 data in U.S. A brief
discussion is presented in Section 6.

##  2 Motivating Data

Our motivating data comes from the COVID tracking project
https://covidtracking.com. State Level COVID-19 Data are recorded for the 50
states plus Washington, DC. For simplicity, we refer to them as “51 states” in
the rest of this paper. Up to June 10th, 2020, United States totally confirmed
2,043,031 cases. 114,533 people died because of COVID-19, and 607,279 people
are recovered from COVID-19. The fatality rate of COVID-19 is $5.6\%$.

Figure 1: Exploratory Analysis of COVID-19 on June 10th

Figure 1 shows state level confirmed numbers, death numbers, incident rate,
and mortality rate. We can see that New York state has the highest confirmed
number, death numbers, and incident rate; Connecticut has the highest
mortality rate among 51 states; Montana has the least confirmed number; Alaska
has the least death number; the incident rate of Hawaii is lowest among 51
states; and Texas has the lowest mortality rate.

##  3 Method

###  3.1 SIRS Model

Compartment epidemic models are widely used to study infectious disease which
spreads through human populations across a large region. SIR model (Kermack
and McKendrick, 1927) has been universally discussed for analyzing the
dynamical evolution of an infectious disease in a large population. SIR model
is extended to SIRS model for imperfect immunity situation (Kermack and
McKendrick, 1932, 1933). For a given time $t$, a fixed population can be split
into three compartments: $S(t)$, $I(t)$, and $R(t)$, which denotes the number
of susceptibles, the number of infectious, and the number of “recovereds”
(which includes deaths), respectively. The dynamical process of SIRS model can
be written as following nonlinear ordinary differential equations of three
given compartments

$\begin{split}&\frac{dS}{dt}=-\beta SI/N+\phi R,\\\ &\frac{dI}{dt}=\beta
SI/N-\gamma I,\\\ &\frac{dR}{dt}=\gamma I-\phi R,\end{split}$ (1)

where $\beta$ denotes the average rate of contact per unit time multiplied by
the probability of disease transmission per contact between a susceptible and
an infectious subject, $\gamma$ denotes the rate of “recovery” per unit time,
which is the rate at which infectious individuals are removed from being
infectious due to recovery (or death), and $\phi$ denotes the rate of loss of
immunity of recovered individuals per unit time, which is the rate at which
recovered individuals become susceptible again (Anderson et al., 1992; Zhuang
and Cressie, 2014). By adding the equations in (1), we notice that

$\frac{dS}{dt}+\frac{dI}{dt}+\frac{dR}{dt}=0.$

Thus, the model postulates a fixed total population without entry and exits of
demographic type. For example, there are no births or deaths caused by other
than the disease we study in a certain time. The SIRS model assume the sum of
all three compartments together is constant within a short period of time such
that

$S(t)+I(t)+R(t)=N,$ (2)

where $N$ is a fixed total population. In cases with discrete time
$t=1,\ldots,T$ (in units of days), we have

$\begin{split}S(t+1)=S(t)-\beta S(t)I(t)/N+\phi R(t),\\\ I(t+1)=I(t)+\beta
S(t)I(t)/N-\gamma I(t),\\\ R(t+1)=R(t)+\gamma I(t)-\phi R(t),\end{split}$ (3)

with the same constraints as (2).

Based on the models in (3) and (2) and assmuptions in (Dukic et al., 2012),
the data model of SIRS assumes conditional independent Poisson distributions
evolving at discrete time points. For a given time $t=1,\ldots,T$, the data
models are

$Z_{R}(t)|P_{R}(t)\sim\text{Possion}(N\times P_{R}(t)),$ (4)

and

$Z_{I}(t)|P_{I}(t)\sim\text{Possion}(N\times P_{I}(t)),$ (5)

where $Z_{R}(t)$ and $Z_{I}(t)$ are the observed number of “recovereds”
(includes deaths) and infectious individuals at time $t$, respectively; $N$ is
known total number of population and $Z_{S}(t)=N-Z_{I}(t)-Z_{R}(t)$;
$P_{R}(t)$ and $P_{I}(t)$ are underlying true rates of recovered and
infectious individuals. Thus, our observed data are
$\\{(Z_{R}(t),Z_{I}(t))\\}:t=1,2\ldots,T$. Based on the relationship between
the number of “recovereds”, infectious, and suspects, we have

$P_{R}(t)+P_{I}(t)+P_{S}(t)=1,$ (6)

where $P_{S}(t)$ underlying rate of susceptible individuals.

Similar to (Zhuang and Cressie, 2014), we have following hidden processes:

$\begin{split}P_{R}(t+1)&=P_{R}(t)+\gamma P_{I}(t)-\phi P_{R}(t),\\\
P_{I}(t+1)&=P_{I}(t)+\beta P_{S}(t)P_{I}(t)-\gamma P_{I}(t),\\\
P_{S}(t+1)&=P_{S}(t)-\beta P_{S}(t)P_{I}(t)+\phi P_{R}(t).\end{split}$ (7)

In order to model the hidden uncertainties in SIRS model, we define following
transformation of $P_{R}(t)$, $P_{I}(t)$ and $P_{S}(t)$ based on (6)

$\begin{split}W_{S}(t)&\equiv\log(\frac{P_{S}(t)}{P_{R}(t)}),\\\
W_{I}(t)&\equiv\log(\frac{P_{I}(t)}{P_{R}(t)}).\\\ \end{split}$ (8)

The time-varying process of $W_{R}(t)$ and $W_{I}(t)$ is defined as

$\begin{split}W_{S}(t+1)&=\mu_{S}(t)+\epsilon_{S}(t+1),\\\
W_{I}(t+1)&=\mu_{I}(t)+\epsilon_{I}(t+1),\end{split}$ (9)

where $\epsilon_{S}(t)\sim N(0,\sigma^{2}_{S})$ and $\epsilon_{I}(t)\sim
N(0,\sigma^{2}_{I})$ for $t=1,2,\ldots,T$. Based on (6) and (7), we have

$\begin{split}\mu_{S}(t)=&W_{S}(t)\\\
&+\log(1+\frac{\phi}{\exp(W_{S}(t))}-\frac{\beta\exp(W_{I}(t))}{1+\exp(W_{S}(t))+\exp(W_{I}(t))})\\\
&+\log(\frac{1}{1+\gamma\exp(W_{I}(t))-\phi})\end{split}$ (10)

and

$\begin{split}\mu_{I}(t)=&W_{I}(t)\\\
&+\log(1-\gamma+\frac{\beta\exp(W_{S}(t))}{1+\exp(W_{S}(t))+\exp(W_{I}(t))})\\\
&+\log(\frac{1}{1+\gamma\exp(W_{I}(t))-\phi})\end{split}$ (11)

Based on the transformation in (8), we have our data in (4) and (5) as

$\begin{split}Z_{R}(t)|W_{S}(t),W_{I}(t)\sim\text{Possion}\left(N\times\frac{1}{1+\exp(W_{S}(t))+\exp(W_{I}(t))}\right),\\\
Z_{I}(t)|W_{S}(t),W_{I}(t)\sim\text{Possion}\left(N\times\frac{\exp(W_{I}(t))}{1+\exp(W_{S}(t))+\exp(W_{I}(t))}\right).\end{split}$
(12)

For the simplicity, we refer the model from (9) to (12) as
$\\{(Z_{R}(t),Z_{I}(t),N),t=1,2\ldots,T\\}\sim\text{SIRS}(\beta,\gamma,\phi,\sigma^{2}_{S},\sigma^{2}_{I})$.
Based on the transmission rate and recover rate, the basic reproduction
number, $R_{0}$, can be calculated by

$R_{0}=\frac{\beta}{\gamma}.$ (13)

###  3.2 Heterogeneity Learning

In section 2, our motivating data is at state level in US and we are
interested in whether there are heterogeneity patterns on the parameters of
interest among different states. As an assumption, we believe that different
states might have different parameters, however, some states will share
similar pattern in transmission rate, recovery rate, or loss of immunity rate.
Next, we introduce nonparametric Bayesian methods for heterogeneity learning
of SIRS parameters over $n$ different regions. In this section, we focus on
the the transmission rate $\beta$ for different regions. Recovery rate and
loss of immunity rate can be parameterized in the same way.

Let $z_{1},\ldots,z_{n}\in\\{1,\ldots,k\\}$ denote clustering labels of $n$
regions and $\beta_{1},\ldots,\beta_{n}$ denote the corresponding parameters
in SIRS model for $n$ regions. Our goal is to cluster
$\beta_{1},\ldots,\beta_{n}$ into $k$ clusters with distinct values
$\beta^{*}_{1},\ldots,\beta^{*}_{k}$, which is usually unknown in practice. A
popular solution for unknown $k$ is to introduce the Dirichlet process mixture
prior models (Antoniak, 1974) as following:

$\displaystyle\beta_{i}\sim G,\quad G\sim DP(\alpha G_{0}),$ (14)

where $G_{0}$ is a base measure and $\alpha$ is a concentration parameter. If
a set of values of $\beta_{1},\ldots,\beta_{n}$ are drawn from $G$, a
conditional prior can be obtained by integration (Blackwell et al., 1973):

$\displaystyle
p(\beta_{n+1}\mid\beta_{1},\ldots,\beta_{n})=\dfrac{1}{n+\alpha}\sum_{j=1}^{n}\delta_{\beta_{j}}(\beta_{n+1})+\dfrac{\alpha}{n+\alpha}G_{0}(\beta_{n+1}).$
(15)

Here, $\delta_{\beta_{j}}(\beta_{\ell})=I(\beta_{\ell}=\beta_{j})$ is a point
mass at $\beta_{j}$. We can obtain the following equivalent models by
introducing cluster membership $z_{j}$’s and letting the unknown number of
clusters $k$ go to infinity (Neal, 2000).

$\begin{split}z_{i}\mid\bm{\pi}&\sim\text{Discrete}(\pi_{1},\ldots,\pi_{k}),\\\
\beta^{*}_{c}&\sim G_{0}\\\
\bm{\pi}&\sim\text{Dirichlet}(\alpha/k,\ldots,\alpha/k)\end{split}$ (16)

where $\bm{\pi}=(\pi_{1},\ldots,\pi_{k})$. In addition, the distribution of
$z_{i}$ can be marginally given by a stick-breaking representation
(Sethuraman, 1994) of Dirichlet process (DP) as

$\begin{split}z_{i}&\sim\sum_{h=1}^{\infty}\pi_{h}\delta_{h},\\\
\pi_{h}&=\nu_{h}\prod_{\ell\leq h}(1-\nu_{\ell}),\\\
\nu_{h}&\sim\text{Beta}(1,\alpha),\end{split}$ (17)

where $\delta_{h}$ is the Dirac function with mass at $h$.

However, Miller and Harrison (2018) proved that the DP mixture model produces
extraneous clusters in the posterior leading to inconsistent estimation of the
number of clusters even when the sample size grows to infinity. A modification
of DP mixture model called Mixture of finite mixtures (MFM) model is proposed
to circumvent this issue (Miller and Harrison, 2018):

$\displaystyle k\sim p(\cdot),\quad(\pi_{1},\ldots,\pi_{k})\mid
k\sim\mbox{Dirichlet}(\eta,\ldots,\eta),\quad z_{i}\mid
k,\bm{\pi}\sim\sum_{h=1}^{k}\pi_{h}\delta_{h},\quad i=1,\ldots,n,$ (18)

where $p(\cdot)$ is a proper probability mass function (p.m.f.) on
$\\{1,2,\ldots,\\}$.

Like the stick-breaking representation in (17) of Dirichlet process, the MFM
also has a similar construction. If we choose $k-1\sim\mbox{Poisson}(\lambda)$
and $\eta=1$ in (18), the mixture weights $\pi_{1},\cdots,\pi_{k}$ is
constructed as follows:

  1. 1.

Generate
$\eta_{1},\eta_{2},\cdots\overset{\text{iid}}{\sim}\text{Exp}(\lambda)$,

  2. 2.

$k=\min\\{j:\sum_{i=1}^{j}\eta_{i}\geq 1\\}$,

  3. 3.

$\pi_{i}=\eta_{i}$, for $i=1,\cdots,k-1$,

  4. 4.

$\pi_{k}=1-\sum_{i}^{k-1}\pi_{i}$.

For ease of exposition, we refer the stick-breaking representation of MFM
above as $\text{MFM}(\lambda)$ with default choice of $p(\cdot)$ being
$\mbox{Poisson}(\lambda)$ and $\eta=1$.

###  3.3 Hierarchical Model

In order to allow for simultaneously heterogeneity learning of three
parameters in SIRS model, the MFM prior is introduced for parameters $\beta$,
$\gamma$ and $\phi$ in the SIRS model. Our observed data are
$\\{(Z_{R}(t,\bm{s}_{i}),Z_{I}(t,\bm{s}_{i}),N_{i}):t=1,2,\ldots,T,i=1,2,\ldots,n\\}$,
where $t$ denotes each discrete time point and $i$ denotes each state. The
hierarchical SIRS model with MFM is given as

$\begin{split}&\\{(Z_{R}(t,\bm{s}_{i}),Z_{I}(t,\bm{s}_{i}),N_{i}),t=1,2\ldots,T\\}\sim\text{SIRS}(\beta_{z^{\beta}_{i}},\gamma_{z^{\gamma}_{i}},\phi_{z^{\phi}_{i}},\sigma^{2}_{S,i},\sigma^{2}_{I,i}),i=1,2,\ldots,n\\\
&z^{\beta}_{i}\sim\text{MFM}(\lambda_{\beta}),i=1,2,\ldots,n,\\\
&z^{\gamma}_{i}\sim\text{MFM}(\lambda_{\gamma}),i=1,2,\ldots,n,\\\
&z^{\phi}_{i}\sim\text{MFM}(\lambda_{\phi}),i=1,2,\ldots,n,\\\
&\beta_{z^{\beta}_{i}}\sim G_{\beta},\\\ &\gamma_{z^{\gamma}_{i}}\sim
G_{\gamma},\\\ &\phi_{z^{\phi}_{i}}\sim G_{\phi},\\\
&\sigma^{2}_{S,i},\sigma^{2}_{I,i}\sim\text{IG}(0.01,0.01),i=1,2,\ldots,n,\end{split}$
(19)

where $z^{\beta}$, $z^{\gamma}$, and $z^{\phi}$ denote the cluster assignments
of parameter $\beta$,$\gamma$, and $\phi$, respectively. $G_{\beta}$,
$G_{\gamma}$, and $G_{\phi}$ is the base distribution for parameter
$\beta$,$\gamma$, and $\phi$, respectively. The choices of $G_{\beta}$,
$G_{\gamma}$, and $G_{\phi}$ will be discussed in Section 3.4.

###  3.4 Prior and Posterior

For the hierarchical SIRS model with MFM introduced in Section 3.3, the set of
parameters is denoted as
$\Theta=\\{\beta_{z^{\beta}_{i}},\gamma_{z^{\gamma}_{i}},\phi_{z^{\phi}_{i}},\sigma^{2}_{S,i},\sigma^{2}_{I,i},\lambda_{\beta},\lambda_{\gamma},\lambda_{\phi}:i=1,2\ldots,n\\}$.
To complete the model, we now specify the joint prior distribution for the
parameters. Based on the natural constraints generated by (3), we have
following distribution for bases distribution $G_{\beta}$, $G_{\gamma}$ and
$G_{\phi}$, respectively:

$\begin{split}&\beta_{z_{i}^{\beta}}\sim\text{Uniform}(0,1),\\\
&\gamma_{z^{\gamma}_{i}}\sim\text{Uniform}(0,1),\\\
&\phi_{z^{\phi}_{i}}\sim\text{Uniform}(0,1).\end{split}$ (20)

For the hyperparameters for three MFM processes, we assign gamma prior
$\text{Gamma}(1,1)$ on $\lambda_{\beta},\lambda_{\gamma},\lambda_{\phi}$. With
the joint prior distributions $\pi(\Theta)$, the posterior distribution of
these parameters based on the data
$D=\\{(Z_{R}(t,\bm{s}_{i}),Z_{I}(t,\bm{s}_{i}),N_{i}):t=1,2,\ldots,T,i=1,2,\ldots,n\\}$
is given as

$\pi(\Theta|(Z_{R}(t,\bm{s}_{i}),Z_{I}(t,\bm{s}_{i}),N_{i}):t=1,2,\ldots,T,i=1,2,\ldots,n)\propto
L(D|\Theta)\times\pi(\Theta),$ (21)

where $L(D|\Theta)$ is the full data likelihood given from the model (9) to
(12). The analytical form of the posterior distribution of
$\pi(\Theta|(Z_{R}(t,\bm{s}_{i}),Z_{I}(t,\bm{s}_{i}),N_{i}):t=1,2,\ldots,T,i=1,2,\ldots,n)$
is unavailable. Therefore, we carry out the posterior inference using the MCMC
sampling algorithm to sample from the posterior distribution and then obtain
the posterior estimates of the unknown parameters. Computation is facilitated
by the nimble package (de Valpine et al., 2017) in R which generates C++ code
for faster computation.

###  3.5 Group Inference via MCMC Samples

After obtaining posterior samples, an important task is do inference on
posterior samples. Using posterior mean or posterior median for grouping label
$\bm{z}$ is not suitable. Instead, inference on the clustering configurations
is obtained employing the modal clustering method of (Dahl, 2006). The
inference is based on the membership matrices of posterior samples,
$B^{(1)},\ldots,B^{(M)}$, where $B^{(t)}$ for the $t$-th post-burn-in MCMC
iteration is defined as:

$\displaystyle
B^{(t)}=[B^{(t)}(i,j)]_{i,j\in\\{1:n\\}}=1(z_{i}^{(t)}=z_{j}^{(t)})_{n\times
n},~{}~{}t=1,\ldots,M,$ (22)

Here $1(\cdot)$ denotes the indicator function, which means $B^{(t)}(i,j)=1$
indicates observations $i$ and $j$ are in the same cluster in the $t$-th
posterior sample post burn-in. After obtaining the membership matrices of the
posterior samples, a Euclidean mean for membership matrices is calculated by:

$\overline{B}=\frac{1}{M}\sum_{t=1}^{M}B^{(t)}.$

Based on $\overline{B}$ and $B^{(1)},\ldots,B^{(M)}$, we find the iteration
with the least squares distance to $\overline{B}$ as

$\displaystyle
C_{L}=\text{argmin}_{t\in(1:M)}\sum_{i=1}^{n}\sum_{j=1}^{n}\\{B(i,j)^{(t)}-\overline{B}(i,j)\\}^{2}.$
(23)

The estimated parameters, together with the cluster assignments $\bm{z}$, are
then extracted from the $C_{L}$-th post burn-in iteration.

##  4 Simulation

In this section, we investigate the performance of the hierarchical SIRS model
with MFM from a variety of measures.

###  4.1 Simulation Settings and Evaluation Metrics

In order to mimic the real dataset we analyze, we choose $n=51$ and the
population for each location is assigned as the real data population for $51$
states. The time length $T$ equals 30 for all the simulation replicates. The
total number of replicates in our simulation study is 100. For each parameter,
we have two different groups and we set the true values of the parameters
$\beta_{1}=0.06,\beta_{2}=0.6$, $\phi_{1}=0.06,\phi_{2}=0.6,$, and
$\gamma_{1}=0.06,\gamma_{2}=0.6$. We randomly assign the labels to $51$
locations and fix them over 100 replicates. The grouping labels for three
parameters is given in Figure 2.

Figure 2: Grouping Labels for $\beta$, $\gamma$, and $\phi$

For each replicates, we have $15,000$ iterations MCMC samples and have first
$5,000$ iterations burn-in in order to obtain samples from every 5th iteration
of the last $10,000$ iterations.

The performance of the posterior estimates of parameters were evaluated by the
mean bias (MB) and the mean standard deviation (MSD) in the following ways,
take $\beta$ as an example:

$\displaystyle\text{MAB}=\frac{1}{100}\sum_{r=1}^{100}\bigg{\\{}\frac{1}{n}\sum_{i=1}^{n}\widehat{\beta}^{r}(\bm{s}_{i})-\beta(\bm{s}_{i})\bigg{\\}},$
$\displaystyle\text{MSD}=\sqrt{\frac{1}{100}\sum_{r=1}^{100}\bigg{\\{}\frac{1}{n}\sum_{i=1}^{n}\left(\widehat{\beta}^{r}(\bm{s}_{i})-\overline{\widehat{\beta}}(\bm{s}_{i})\right)^{2}}\bigg{\\}},$

where $\overline{\widehat{\beta}}(\bm{s}_{i})$ is the mean of the posterior
estimate over 100 replicates.

For clustering estimation evaluation, the estimated number of clusters
$\widehat{K}$ for each replicate is summarized from the MCMC iteration picked
by Dahl’s method. Rand Index (RI; Rand, 1971) is applied to evaluate cluster
configuration. The RI is calculated by R-package fossil (Vavrek, 2011). A
higher value of the RI represents higher accuracy of clustering The average RI
(MRI) was calculated as the mean of RIs over the 100 replicates.

###  4.2 Simulation Results

The parameter estimation performance and clustering performance results are
shown in Table 1 and Table 2.

Table 1: Estimation Performance for Simulation Study  Parameter | MB | MSD
---|---|---
$\beta_{1}$ | 0.008 | 0.021
$\beta_{2}$ | -0.072 | 0.152
$\gamma_{1}$ | 0.007 | 0.017
$\gamma_{2}$ | -0.068 | 0.151
$\phi_{1}$ | 0.012 | 0.023
$\phi_{2}$ | -0.069 | 0.149
Table 2: Grouping Performance for Simulation Study  Parameter | MRI | S.D of RI | $\widehat{K}$ | S.D. of $\widehat{K}$
---|---|---|---|---
$\beta$ | 0.854 | 0.058 | 2.12 | 0.33
$\gamma$ | 0.857 | 0.057 | 2.33 | 0.55
$\phi$ | 0.847 | 0.059 | 2.31 | 0.54

From the results shown in Table 1, we see that the MABs and MSDs of the
parameters are both within a reasonable range. In general, performance of
posterior estimates achieve a good target.

And we see that our proposed methods successfully recover the number of groups
and grouping labels within a reasonable range for all three parameters from
Table 2. Average rand index for all parameters around 0.85 indicate our
proposed method truly recovers the group labels for all three parameters. The
mean of the estimated number of groups is close to true number of groups over
100 replicates.

##  5 Real Data Analysis

###  5.1 30-Day Analysis from April 1st

We analyze 30-Day data from April 1st, 2020. The reason why we analyze this
time period data is that U.S. Government announced the suspension of entry as
immigrants and nonimmigrants of certain additional persons who pose a risk of
transmitting corona-virus https://www.whitehouse.gov/presidential-actions/ on
March 11th, 2020. From the April 1st, we can assume that there are very
limited imported cases from outside U.S.. We analyze 30-day data based on the
model in (19) and use the priors discussed in Section 3.4. We run $50,000$
MCMC iterations and burnin the first $20,000$ iterations in order to obtain
samples from every 10th iteration of the last $30,000$ iterations. The group
labels are obtained by Dalh’s method in Section 3.5.

For $\beta$, one group is identified. $\beta=0.079$ with $95\%$ Highest
Probability Density (HPD) interval (Chen and Shao, 1999) $(0.058,0.098)$. For
$\gamma$, three groups are identified with $\gamma_{1}=0.0054$ with HPD
interval $(0.0021,0.0207)$, $\gamma_{2}=0.0419$ with HPD interval
$(0.0022,0.0609)$ and $\gamma_{3}=0.0164$ with HPD interval$(0.0035,0.0241)$.
Thirty three states including Alabama, Arizona, California, Colorado,
Connecticut, District of Columbia, Florida, Georgia, Idaho, Illinois, Indiana,
Kansas, Louisiana, Massachusetts, Michigan, Mississippi, Missouri, Nebraska,
Nevada, New Jersey, North Carolina, Ohio, Oregon, Pennsylvania, Rhode Island,
South Carolina, Texas, Utah, Vermont, Virginia, Washington, West Virginia,
Wisconsin, belong to group one. Eleven states including Arkansas, Hawaii,
Iowa, Maine, Minnesota, New Hampshire, North Dakota, Oklahoma, South Dakota,
Tennessee, Wyoming , belong to group 2. And seven states including Alaska,
Delaware, Kentucky, Maryland, Montana, New Mexico, New York, belong to group
3. The grouping labels for $\gamma$s are shown in Figure 3.

For $\phi$, one group is identified. $\phi=0.0015$ with HPD interval
$(1.181\times 10^{-7},0.0047)$.

Figure 3: Group Labels for $\gamma$s of April 1st data

With the estimated values of $\beta$ and $\gamma$, the basic reproduction
number, $R_{0}$, is calculated among different states. The values of $R_{0}$
among different states are shown in Figure 4.

Figure 4: $R_{0}$ for 51 States from April 1st

###  5.2 30-Day Analysis from May 1st

The second time period we analyze is from May 1st, 2020. Other settings are
same with previous analysis.

For $\beta$, one group is identified. $\beta=0.0042$ with $95\%$ Highest
Probability Density (HPD) interval (Chen and Shao, 1999) $(3.056\times
10^{-8},0.1083)$. Compared with previous 30-day data, in this time period, the
transmission rate decrese a lot. For $\gamma$, two groups are identified with
$\gamma_{1}=0.0381$ with HPD interval $(0.0048,0.3713)$ and
$\gamma_{2}=0.0007$ with HPD interval $(0.0003,0.0013)$. Two states, Oregon
and Vermont, are identified in group 1. Other states are identified in group
2. For $\phi$, one groups is identified. $\phi=0.0006$ with HPD interval
$(2.747\times 10^{-7},0.0026)$.

With the estimated values of $\beta$ and $\gamma$, the basic reproduction
number, $R_{0}$, is calculated among different states. There are two different
groups for the basic reproduction number. One group include Oregon and Vermont
with $R_{0}=0.1102$. The other group includes other 49 states with
$R_{0}=5.4619$. Comparing to the 30 days period from April 1st, we can see a
decrease for $R_{0}$ in general.

##  6 Discussion

In this paper, we develop a nonparametric Bayesian heterogeneity learning
method for SIRS model based on Mixture of Finite Mixtures model. This
statistical framework was motivated by the heterogeneity of COVID-19 pattern
among different regions. Our simulation results indicate that the proposed
method can recover the heterogeneity pattern of parameters among different
regions. Illustrated by the analysis of COVID-19 data in U.S., our proposed
methods reveal the heterogeneity pattern among different states.

In addition, three topics beyond the scope of this paper are worth further
investigation. First, we can add spatially dependent structure (Hu et al.,
2020; Zhao et al., 2020) on the heterogeneity of different states. Second, our
model assumes parameters are constant over a certain time period. Building
heterogeneity learning model with time varying coefficients is an interesting
future work. Finally, proposing a measurement error correction for SIRS
devotes another interesting future work.

## References

  * Anderson et al. (1992) Anderson, R. M., B. Anderson, and R. M. May (1992).  Infectious diseases of humans: dynamics and control.  Oxford university press.
  * Antoniak (1974) Antoniak, C. E. (1974).  Mixtures of Dirichlet processes with applications to Bayesian nonparametric problems.  The Annals of Statistics 2(6), 1152–1174.
  * Blackwell et al. (1973) Blackwell, D., J. B. MacQueen, et al. (1973).  Ferguson distributions via Pólya urn schemes.  The Annals of Statistics 1(2), 353–355.
  * Chen and Shao (1999) Chen, M.-H. and Q.-M. Shao (1999).  Monte carlo estimation of bayesian credible and hpd intervals.  Journal of Computational and Graphical Statistics 8(1), 69–92.
  * Chen et al. (2020) Chen, Y.-C., P.-E. Lu, C.-S. Chang, and T.-H. Liu (2020).  A time-dependent sir model for covid-19 with undetectable infected persons.  arXiv preprint arXiv:2003.00122.
  * Dahl (2006) Dahl, D. B. (2006).  Model-based clustering for expression data via a dirichlet process mixture model.  Bayesian inference for gene expression and proteomics 4, 201–218.
  * de Valpine et al. (2017) de Valpine, P., D. Turek, C. J. Paciorek, C. Anderson-Bergman, D. T. Lang, and R. Bodik (2017).  Programming with models: writing statistical algorithms for general model structures with nimble.  Journal of Computational and Graphical Statistics 26(2), 403–413.
  * Dukic et al. (2012) Dukic, V., H. F. Lopes, and N. G. Polson (2012).  Tracking epidemics with google flu trends data and a state-space seir model.  Journal of the American Statistical Association 107(500), 1410–1426.
  * Geng et al. (2019) Geng, J., A. Bhattacharya, and D. Pati (2019).  Probabilistic community detection with unknown number of communities.  Journal of the American Statistical Association 114(526), 893–905.
  * Geng et al. (2019) Geng, J., W. Shi, and G. Hu (2019).  Bayesian nonparametric nonhomogeneous poisson process with applications to usgs earthquake data.  arXiv preprint arXiv:1907.03186.
  * Hethcote (2000) Hethcote, H. W. (2000).  The mathematics of infectious diseases.  SIAM review 42(4), 599–653.
  * Hu et al. (2020) Hu, G., J. Geng, Y. Xue, and H. Sang (2020).  Bayesian spatial homogeneity pursuit of functional data: an application to the us income distribution.  arXiv preprint arXiv:2002.06663.
  * Jo et al. (2020) Jo, H., H. Son, S. Y. Jung, and H. J. Hwang (2020).  Analysis of covid-19 spread in south korea using the sir model with time-dependent parameters and deep learning.  medRxiv.
  * Kermack and McKendrick (1927) Kermack, W. O. and A. G. McKendrick (1927).  A contribution to the mathematical theory of epidemics.  Proceedings of the royal society of london. Series A, Containing papers of a mathematical and physical character 115(772), 700–721.
  * Kermack and McKendrick (1932) Kermack, W. O. and A. G. McKendrick (1932).  Contributions to the mathematical theory of epidemics. ii.—the problem of endemicity.  Proceedings of the Royal Society of London. Series A, containing papers of a mathematical and physical character 138(834), 55–83.
  * Kermack and McKendrick (1933) Kermack, W. O. and A. G. McKendrick (1933).  Contributions to the mathematical theory of epidemics. iii.—further studies of the problem of endemicity.  Proceedings of the Royal Society of London. Series A, Containing Papers of a Mathematical and Physical Character 141(843), 94–122.
  * Miller and Harrison (2013) Miller, J. W. and M. T. Harrison (2013).  A simple example of Dirichlet process mixture inconsistency for the number of components.  In Advances in Neural Information Processing Systems, pp. 199–206.
  * Miller and Harrison (2018) Miller, J. W. and M. T. Harrison (2018).  Mixture models with a prior on the number of components.  Journal of the American Statistical Association 113(521), 340–356.
  * Neal (2000) Neal, R. M. (2000).  Markov chain sampling methods for Dirichlet process mixture models.  Journal of Computational and Graphical Statistics 9(2), 249–265.
  * Rand (1971) Rand, W. M. (1971).  Objective criteria for the evaluation of clustering methods.  Journal of the American Statistical association 66(336), 846–850.
  * Read et al. (2020) Read, J. M., J. R. Bridgen, D. A. Cummings, A. Ho, and C. P. Jewell (2020).  Novel coronavirus 2019-ncov: early estimation of epidemiological parameters and epidemic predictions.  MedRxiv.
  * Sethuraman (1994) Sethuraman, J. (1994).  A constructive definition of Dirichlet priors.  Statistica Sinica, 639–650.
  * Sun et al. (2020) Sun, H., Y. Qiu, H. Yan, Y. Huang, Y. Zhu, J. Gu, and S. X. Chen (2020).  [discussion paper] tracking reproductivity of covid-19 epidemic in china with varying coefficient sir model.  Journal of Data Science, 2.
  * Tang et al. (2020) Tang, B., X. Wang, Q. Li, N. L. Bragazzi, S. Tang, Y. Xiao, and J. Wu (2020).  Estimation of the transmission risk of the 2019-ncov and its implication for public health interventions.  Journal of clinical medicine 9(2), 462.
  * Vavrek (2011) Vavrek, M. J. (2011).  Fossil: palaeoecological and palaeogeographical analysis tools.  Palaeontologia Electronica 14(1), 16.
  * Wu et al. (2020) Wu, J. T., K. Leung, and G. M. Leung (2020).  Nowcasting and forecasting the potential domestic and international spread of the 2019-ncov outbreak originating in wuhan, china: a modelling study.  The Lancet 395(10225), 689–697.
  * Zhao et al. (2020) Zhao, P., H.-C. Yang, D. K. Dey, and G. Hu (2020).  Bayesian spatial homogeneity pursuit regression for count value data.  arXiv preprint arXiv:2002.06678.
  * Zhuang and Cressie (2014) Zhuang, L. and N. Cressie (2014).  Bayesian hierarchical statistical sirs models.  Statistical Methods & Applications 23(4), 601–646.