\section{Introduction and preliminaries}
We consider the so-called \emph{simple games} and
the computational complexity of representing them as unions or intersections of weighted majority games.
Simple games and its dimension, as well as, weighted majority games, were defined by Taylor and Zwicker~\cite{TaZw99}.
Later, Freixas and Marciniak~\cite{FM10} introduced a new concept, the codimension of simple games.

A {\em simple game} is a tuple $\Gamma=(N,\cW)$, where $N$ is a finite set of {\em players} and $\cW\subseteq\cP(N)$ is a monotonic family of subsets of $N$.
Furthermore, \emph{its dual} $\Gamma^*=(N,\cW^*)$ is the game such that
$\cW^*=\{S\subseteq{}N\,:\,N\setminus{}S\not\in{}\cW\}$.
$\Gamma$ is said to be \emph{self-dual} if $\Gamma=\Gamma^*$.
Note that $(\Gamma^*)^*=\Gamma$.
Given two simple games $\Gamma_1=(N_1,\cW_1)$ and $\Gamma_2=(N_2,\cW_2)$, they are \emph{equivalent}
if $N_1=N_2$ and $\cW_1=\cW_2$.
The subsets of $N$ are called {\em coalitions}, the set $N$ is the {\em grand coalition} and each $X\in\cW$ is a {\em winning coalition}.
The complement of the family of winning coalitions is the family of {\em losing coalitions} $\cL$, i.e., $\cL=\cP(N)\setminus\cW$.
Any of those set families determine uniquely the game $\Gamma$ and constitute one of the usual forms of representation for simple games~\cite{TaZw99}, although the size of the representation is not, in general, polynomial in the number of players~\cite{MRS15b}.

A simple game $\Gamma$ is a {\em weighted majority game} (\WMG) if it admits a representation by means of  $n+1$ nonnegative real numbers $[q;w_1,\ldots,w_n]$ such that $S\in\cW
\iff w(S)\ge{}q$ where, for each coalition $S\subseteq{}N$, $w(S)=\sum_{i\in{}S}w_i$. The number $q$ is called the quota and $w_i$ the weight of the player $i$. It is well known that any \WMG admits a representation with integer numbers.

The {\em dimension} of a simple game $\Gamma$ is the least $k$ such that there exists \WMG{}s $\Gamma_1,\ldots,\Gamma_k$ such that $\Gamma=\Gamma_1\cap\ldots\cap\Gamma_k$.
On the other hand, the {\em codimension} of a simple game $\Gamma$ is the least $k$ such that there exists \WMG{}s $\Gamma_1,\ldots,\Gamma_k$ such that $\Gamma=\Gamma_1\cup\ldots\cup\Gamma_k$.

Many theoretical results and examples about dimension and codimension~\cite{OlKuMo16,KuNa15,FM10,FrPu08,FP01,TaZw99,FrPu98}
have been constantly appearing during the last years,
as well as, computational complexity results~\cite{MOS16,FMOS11,DeWo06}.
We present some results that will be used in Section~\ref{sec:comcom}.
\begin{lemma}
  \label{lemma_upper_bound_maximal_losing}
  The dimension of a simple game $v$ is bounded above by $\left|\mathcal{L}^M\right|$ and the co-dimension is bounded above by
  $\left|\cW^m\right|$.
\end{lemma}
\remove{
\begin{proof}
  For each coalition $S\in \mathcal{L}^M$ we set $q^S=1$, $w_i^S=0$ for all $i\in S$ and $w_i^S=1$ otherwise. Note that
  $\mathcal{L}^M\neq \emptyset$ since $\emptyset$ is a losing coalition. With this $w^S(S)=0<q^S$. However, for all $T\subseteq N$
  with $T\not\subseteq S$ we have $w(T)\ge 1=q^S$. Thus, we have $v=\bigwedge_{S\in\mathcal{L}^M} \left[q^S;w_1^S,\dots,w_n^S\right]$.
  Similarly, for each $S\in\cW^n$ we set $\widetilde{q}^S=|S|$, $\widetilde{w}_i^S=1$ for all $i\in S$ and $w_i^S=0$ otherwise.
  Note that $\cW^m\neq \emptyset$ since $N$ is a winning coalition. With this $\widetilde{w}^S(S)=\widetilde{q}^S$. However, for all
  $T\subseteq N$ with $S\not\subseteq T$ we have $w(T)<\widetilde{q}^S$. Thus, we have
  $v=\bigvee_{S\in\cW^m} \left[\widetilde{q}^S;\widetilde{w}_1^S,\dots,\widetilde{w}_n^S\right]$.
\end{proof}
}

\begin{lemma}
\label{lem:game-dual}
Let $\Gamma$ be a simple game.
 $\Gamma$ is the intersection of $t$ weighted games if and only if
 $\Gamma^*$ is the union of $t$ weighted games.
Furthermore a representation, as union (intersection), of $\Gamma^*$ can be obtained from a representation, as intersection (union), of $\Gamma$  in polynomial time.
 Moreover, dim($\Gamma$)=codim($\Gamma^*$),
 and if $\Gamma$ is self-dual then $dim(\Gamma)=codim(\Gamma)$.
\end{lemma}

Note that the converse statement of the last sentence is not true in general as there are weighted games which are not self-dual.

\section{Computational complexity of related problems}
\label{sec:comcom}

First, we present a simple game with $2n$ players, dimension $n$ and codimension $2^{n-1}$.
Other examples of simple games with $2n$ players and dimension $2^{n-1}$ can be found in~\cite{TaZw99,OlKuMo16}.

\begin{example}
\label{ex:dim-codim}
Given a positive integer $n$,
Freixas and Marcinicak (Theorem 2 of~\cite{FM10})
define a simple game
with $2n$ players and dimension $n$.
Let $\Gamma=(N,\cW)$ be a simple game
defined by $N=\{1,2,\ldots,2n\}$ and $S\in{}\cW$ iff
$S\cap\{2i-1,2i\}\neq\emptyset$, $i\in\{1,\ldots,n\}$,
then $\Gamma$ has dimension $n$,
\remove{
\[
\begin{array}{rl}
 \Gamma = & \left[1; 1,1,0,0,\ldots,0,0\right] \cap \left[1; 0,0,1,1,0,0,\ldots,0,0\right] \cap \ldots \cap\\[.25cm]
          & \left[1;0,0,\ldots,0,0,1,1,0,0\right] \cap \left[1; 0,0,\ldots,0,0,1,1\right]
\end{array}
\]
}
\[
 \Gamma = \left[1; 1,1,0,0,\ldots,0,0\right] \cap \left[1; 0,0,1,1,0,0,\ldots,0,0\right]
          \cap \ldots \cap \left[1; 0,0,\ldots,0,0,1,1\right]
\]
As $S$ is a winning coalition in $\Gamma^*$ iff
$N\setminus{}S$ is a losing coalition in $\Gamma$,
$ \Gamma^*=(N, \cW_1^* \cup \ldots \cup \cW_n^*)$,
where $\cW_i^*=\{S\subseteq{}N\,:\, \{2i-1,2i\}\subseteq{}S\}$, $i\in{}\{1,\ldots,n\}$.
As $\Gamma^*$ is a composition of $n$ unanimity games,
$\Gamma^*$ has dimension $2^{n-1}$~\cite{FP01} and
$\Gamma$ has codimension $2^{n-1}$ (by Lemma~\ref{lem:game-dual}).
\end{example}

\begin{proposition}
Given a simple game $\Gamma$ as union (intersection) of weighted games,
computing a representation of $\Gamma$ as intersection (union) of weighted games requires exponential time.
\end{proposition}

The complexity of
several problems about representations of simple games as intersections of \WMG{}s were analyzed in~\cite{DeWo06}. We provide here a new reduction from the  NP-hard \emph{Subset Sum Problem (\SSP)}. Our reduction differs in the fact that  for the associated game  we know both the dimension and the codimension.  Recall that in the \SSP on input
 a finite set $A=\{a_1,\ldots,a_n\}$ of positive integers and a positive integer $b$ we want to know whether  there is a subset $A'\subseteq{}A$ such that $\sum_{a_i\in{}A'}a_i=b$.

Let  $I=(b; a_1, a_2, \ldots, a_n)$ be an instance of the \SSP, $d>1$,
and let $\Gamma(I,d)$ be the  game defined on $n+d$ players by the intersection of the $d$ \WMG{}s:
\[
 \begin{array}{l}
  [3b+1; 3 a_1, \ldots, 3 a_n, 1, 1, 0, \stackbin{2(d-1)}{\ldots}, 0],\
  [3b+1; 3 a_1, \ldots, 3 a_n, 0, 0, 1, 1, 0, \stackbin{2(d-2)}{\ldots}, 0],\\
  \ldots,\  [3b+1; 3 a_1, \ldots, 3 a_n, 0, \stackbin{2(d-1)}{\ldots}, 0, 1, 1].
 \end{array}
\]

\begin{lemma}
\label{prop:inter}\label{lem:red}
Let $d>1$. When  $I$ is a yes instance of \SSP then $dim(\Gamma(I,d))=d$ and $codim(\Gamma(I,d))=2^{d}$, otherwise,   $dim(\Gamma(I,d))=codim(\Gamma(I,d))=1$.
\end{lemma}
\begin{proof}
Let $X\subset\{1,\ldots,n\}$,  $E=\{n+1,\dots, n+2d\}$ and $Y\subset E$.  Let $\Gamma$ be the game given in Example 1 taking $E$ as the set of players.
Observe that if $\sum_{i\in{}X}a_i> b$,  $\sum_{i\in{}X}a_i\ge{}b+1$, thus  $X \cup Y$ is a winning coalition in  $\Gamma(I,d)$.
When $\sum_{i\in{}X}a_i< b$, $\sum_{i\in{}X}a_i \le{}b-1$, thus  $X \cup Y$ is a losing coalition in  $\Gamma(I,d)$.
In the case that  $\sum_{i\in{}X}a_i= b$,   $X \cup Y$ is a winning coalition in  $\Gamma(I,d)$ iff $Y$ is a winning coalition in $\Gamma$.

When  $I$ is a yes instance of \SSP it follows that  $dim(\Gamma(I,d)) =d$ as otherwise $\Gamma$ will have a smaller dimension. On the other hand, it is easy to see that $\Gamma(I,d)^*$ is the composition of $d+1$ unanimity games and therefore $dim(\Gamma(I,d)^*) =2^d$ according to  \cite{FP01}, thus $codim(\Gamma(I,d)) =2^d$.
When $I$ is a  no instance, there is  no $X\subseteq \{1,\dots,n\}$  for which  $\sum_{i\in{}X}a_i= b$. Therefore  $X\cup Y$ is winning in $\Gamma(I,d)$ iff $\sum_{i\in{}X}a_i\geq  b+1$.
Thus $\Gamma(I,d)$ is the game $[b+1; a_1,\dots,a_n, 0, \stackbin{2d}{\ldots}, 0]$ and we conclude that $dim(\Gamma(I,d))=codim(\Gamma(I,d))=1$.
\end{proof}

Combining lemmas~\ref{lem:red} and \ref{lem:game-dual}   we can prove  the following results.

\begin{proposition}
\label{prop:inter}
Let $d_1$ and $d_2$ be two integers with $1\le{}d_2<d_1$.
Then the problem of deciding whether the union of $d_1$ given \WMG{}s can also be represented as the union of $d_2$ \WMG{}s is NP-hard.
\end{proposition}

\begin{proposition}
\label{prop:codim1}
Let $d_1$ and $d_2$ be two integers with $1\le{}d_1,d_2$.
Then the problem of deciding whether the intersection (union) of $d_1$ given \WMG{}s can also be represented
as the union (intersection) of $d_2$ \WMG{}s is NP-hard.
\end{proposition}

As a consequence of the previous results, given a simple game $\Gamma$ as union or intersection of \WMG{}s,
to compute $dim(\Gamma)$, $codim(\Gamma)$ or deciding whether $\Gamma$ is weighted are  NP-hard problems.

Recall that two game representations are said to be \emph{equivalent} whenever  the represented games have the same set of winning coalitions.
We can extend several results on equivalence problems from ~\cite{ElkindGGW08} to games given as unions of \WMG, in particular we have.

\begin{proposition}
\label{theo:Elkind4}
Checking whether a given  union of \WMG{}s
is equivalent to a given  union of \WMG{}s is co-NP-complete,
even if all weights are equal to 0 or 1.
\end{proposition}

\section{Future work}
It remains open to exhaustively classify the dimension and codimension of all complete simple game  up to $n$ players.
Some bounds about dimension are given by Freixas and Puente~\cite{FP01} and Olsen \emph{et al.}~\cite{OlKuMo16}.
As well as to find complete simple games with \emph{small} dimension (codimension), but
with \emph{large} codimension (dimension), and to construct analytical examples with specific dimension and codimension.
It is also interesting to find real simple games with \emph{large} dimension or codimension as the example given by Kurz and Napel~\cite{KuNa15}.

\section*{References}


\section{Introduction}\label{section: introduction}

Many popular measures of pairwise dependence, for example Kendall's
tau \citep{Kendall38} and Spearman's rho \citep{Spearman04}, have the
undesirable property that they may be zero even when the two
considered random variables $X$ and $Y$ are dependent.  Addressing
this weakness, \cite{BergsmaDassios14} introduced a new rank-based
correlation measure $\tau^*$, which, under mild conditions on the
joint distribution of $(X,Y)$, is zero if and only if $X$ and $Y$ are
independent.  Where Kendall's tau is defined in terms of concordance
and discordance of two independent copies of $(X,Y)$, the new $\tau^*$
is based on similar notions of concordance and discordance for four
independent copies of $(X,Y)$.  While a na\"{\i}ve computation of
$t^*$, the empirical version of $\tau^*$, thus requires $O(n^4)$ time
for a sample of size $n$, it was recently shown that this
computational burden can be reduced to $O(n^2\log(n))$
\citep{WeihsEtAl15}.  As $t^*$ is now computable for larger sample
sizes, understanding its asymptotic behavior becomes a problem of
practical interest and has the potential to yield simple tests of
independence that avoid Monte Carlo approximation of p-values.

We introduce the statistic $t^*$ in Section~\ref{section:
  preliminaries}, where we also review background on U-statistics.  In
Section~\ref{sec:degeneracy}, we clarify that $t^*$ is a degenerate
U-statistic under the null hypothesis that the sample is generated
under independence.  We also prove that in certain settings degeneracy
occurs only under independence.  In Section \ref{section: asymptotics
  under the null}, we use the asymptotic theory of degenerate
U-statistics to derive an explicit representation of the asymptotic
distribution of $t^*$ when the sample is generated under independence
and with marginals that are continuous or discrete.  The asymptotic
distribution takes the form of a Gaussian chaos; specifically, we find
a (in some cases infinite) sum of scaled and centered chi-square
distributions.  Simulations in Section \ref{section: simulations} then
demonstrate how the large-sample theory can be leveraged to perform
tests of independence and compute power.  Indeed, asymptotic
distributions are found to give accurate approximations for sample
sizes as small as $n=80$.  We end with a discussion in Section \ref{section: discussion}.

\section{Preliminaries}\label{section: preliminaries}

\subsection{The $t^*$ statistic} \label{subsection: t*} Let
$(x_1,y_1),\ldots,(x_n,y_n)$ be a sample of points in $\mathbb{R}^2$.
The empirical version of the Bergsma-Dassios sign covariance is
the statistic
\begin{align}\label{equation: test statistics}
t^* := \frac{(n-4)!}{n!} \sum_{\substack{1\leq i,j,k,l \leq n \\
    i,j,k,l~\text{distinct}}}
a(x_i,x_j,x_k,x_l) a(y_i,y_j,y_k,y_l),
\end{align}
where
\begin{multline}
  \label{eq:def-a}
a(z_1,z_2,z_3,z_4) =\\ I(z_1,z_3<z_2,z_4) + I(z_1,z_3>z_2,z_4) -
I(z_1,z_2<z_3,z_4) - I(z_1,z_2>z_3,z_4).
\end{multline}
Here we use $I(\cdot)$ to denote the indicator function and $a,b < c,d$ is shorthand for $\max(a,b) < \min(c,d)$.  As in
\cite{WeihsEtAl15}, we defined $t^*$ in the form of a U-statistic,
whereas \cite{BergsmaDassios14} introduced it as a V-statistic.
Indeed, $t^*$ from~(\ref{equation: test statistics}) is an unbiased
estimator of the sign covariance
\[
\tau^* :=
\Exp\left[a(X_1,X_2,X_3,X_4)a(Y_1,Y_2,Y_3,Y_4)\right]
\]
of \cite{BergsmaDassios14}.  Here, $(X_1,Y_1),...,(X_4,Y_4)$ are
random vectors drawn independently from a given bivariate distribution on
$\mathbb{R}^2$.

\begin{exmp}
  \label{ex:gaussian-taustar}
  Figure~\ref{fig:normal-taustar} shows the values of $\tau^*$ for
  bivariate normal distributions, which we computed by Monte Carlo
  simulation.  The sign covariance $\tau^*$ is an even function of the
  normal correlation $\rho$, and we thus only show values for
  $\rho\in[0,1]$.  For each considered correlation $\rho$ we
  averaged 200 values of $t^*$, each computed from a
  sample of size $n=300$.
\end{exmp}

\begin{figure}[t]
  \centering
  \includegraphics[scale=.4]{bivariate_normal_relation}
  \caption{Sign covariance $\tau^*$ of bivariate normal
    distributions.}
  \label{fig:normal-taustar}
\end{figure}

As we explain in the remainder of this subsection, the statistic $t^*$
is based on counting concordant and disconcordant quadruples.

\begin{defn}
Let $(x_1,y_1),...,(x_4,y_4)$ be four points relabelled so that $x_1\leq x_2\leq x_3\leq x_4$.  We say that the points are
\begin{align*}
	\begin{array}{lll}
		\text{\emph{inseparable}} & \mbox{if } \ \
			\begin{array}{@{}l@{}}
				\text{$x_2=x_3$ or there exists a permutation $\pi$ of $\{1,2,3,4\}$} \\
				\text{so that $y_{\pi(1)}\leq y_{\pi(2)}=y_{\pi(3)}\leq y_{\pi(4)}$,}
			\end{array}
	\end{array}
\end{align*}
and if they are not inseparable, then we call them
\begin{align*}
		\begin{array}{lll}
			\text{\emph{concordant}} & \mbox{if } \text{$\max(y_1,y_2) < \min(y_3,y_4)$\ \ or \ $\max(y_3,y_4) < \min(y_1,y_2)$,} \\
			\text{\emph{discordant}} & \mbox{if }  \text{$\max(y_1,y_2) > \min(y_3,y_4)$ and $\max(y_3,y_4) > \min(y_1,y_2) $.}
		\end{array}
\end{align*}
\end{defn}
The above definitions are mutually exclusive and exhaustive in that
any set of four points in $\mathbb{R}^2$ will be exactly one of inseparable,
concordant, or discordant. Moreover, if the points are drawn from a
bivariate distribution with continuous marginals then they will be
almost surely concordant or discordant. See Figure 3 of
\cite{BergsmaDassios14} for a visual depiction of concordance and
discordance.

Let $S_4$ be the set of permutations on 4 elements, and for $\pi \in
S_4$ and $(z_1,z_2,z_3,z_4) \in \mathbb{R}^{4}$, let $z_{\pi(1,2,3,4)}
:= (z_{\pi(1)},z_{\pi(2)},z_{\pi(3)},z_{\pi(4)})$.  Introducing the
symmetric function
\begin{align}\label{eq: kernel}
h((x_1,y_1),\ldots,(x_4,y_4)) := \frac{1}{4!} \sum_{\pi \in S_4}
a(x_{\pi(1,2,3,4)}) a(y_{\pi(1,2,3,4)}),
\end{align}
we may rewrite $t^*$ as a sum of permutation invariant terms,
namely,
\begin{align}\label{eq: Bergsma-Dassios U-statistic}
t^* = \frac{1}{\binom n4} \sum_{(i,j,k,l) \in C(n,4)} h\left((x_i,y_i),(x_j,y_j),(x_k,y_k),(x_l,y_l)\right),
\end{align}
where $C(n,4)= \{(i,j,k,l) \,:\,1 \leq i < j < k < l \leq n \}$.  Lemma
1 in \cite{WeihsEtAl15} gives the following result.

\begin{lem} \label{lemma: kernel values}
Let $A = \{(x_1,y_1),(x_2,y_2),(x_3,y_3),(x_4,y_4)\}\subset\mathbb{R}^2$. Then
\begin{align*}
	h((x_1,y_1),..,(x_4,y_4)) = \left\{
		\begin{array}{lll}
			2/3 & \mbox{if } \text{the points in $A$ are concordant,}\\
			-1/3 & \mbox{if }  \text{the points in $A$ are discordant,} \\
			0 & \mbox{if }  \text{the points in $A$ are inseparable.}
		\end{array}
	\right.
\end{align*}
\end{lem}

Equation~(\ref{eq: Bergsma-Dassios U-statistic}) expresses $t^*$ in
the familiar form of a U-statistic with symmetric kernel $h$, and we
proceed to review some of the tools available for the study of
U-statistics.

\subsection{Theory of U-statistics} \label{subsection: theory of
  u-statistics} Let $Z_1,Z_2,\dots$ be i.i.d.\ random variables taking
their values in $\mathbb{R}^d$ with $d\geq 1$.  Let $k:\left(\mathbb{R}^d\right)^m
\to \mathbb{R}$ be a \emph{kernel function} invariant to permutation of its
$m$
arguments. For $n\geq m$, the \emph{U-statistic with kernel $k$} is
the statistic
\begin{align}
  \label{eq:Un}
	U_n := \frac{1}{{n \choose m}} \sum_{(i_1,...,i_m)\in C(n,m)} k(Z_{i_1},...,Z_{i_m}),
\end{align}
where $C(n,m) = \{(i_1,...,i_m)\in \{1,...,n\}^m: i_1 < i_2 < ... <
i_m\}$.  Note that $\Exp[U_n] = \Exp[k(Z_1,...,Z_m)]$ so that $U_n$ is
an unbiased estimator of $\theta:=\Exp[k(Z_1,...,Z_m)]$.

Of central importance in determining the asymptotics of U-statistics
are the functions
\begin{equation}
  \label{eq:k_i}
  k_i(z_1,...,z_i) =  \Exp[k(z_1,...,z_i, Z_{i+1},\dots,Z_{m})], \quad i
=1,\dots,m,
\end{equation}
 and their variances
\begin{equation}
  \label{eq:sigma_i^2}
  \sigma_i^2 = \Var[k_i(Z_1,...,Z_i)], \quad i=1,...,m.
\end{equation}
It is well known that $\sigma_1^2 \leq \sigma_2^2 \leq ...\leq
\sigma_m^2$.  In particular, if $\sigma_m^2$ is finite then so are all
other $\sigma_i^2$.  We now recall two theorems on the large-sample
distribution of the U-statistic $U_n$ \cite[Chapter 5]{Serfling:1980}.

\begin{thm} \label{theorem: normal asymptotics}

  If the kernel $k$ of the statistic $U_n$ from~(\ref{eq:Un})
  has variance $\sigma_m^2 < \infty$, then
\begin{align*}
	\sqrt{n}(U_n-\theta) \overset{d}{\to} N(0, m^2 \sigma_1^2).
\end{align*}
\end{thm}

If $\sigma_1^2 = 0$, then the Gaussian limit is degenerate, and we
have $\sqrt{n}(U_n-\theta) \overset{p}{\to} 0$.  Indeed, if
$\sigma_1^2 = 0$ and $\sigma_2^2 >0$, then scaling $U_n$ by a factor
of $n$ results in a non-Gaussian asymptotic distribution.  To present
this result, we write $\chi^2_1$ for the chi-square distribution with
one degree of freedom and define $A_k$ to be the operator that acts
via $g(\cdot) \mapsto \Exp[(k_2(\cdot, Z_1) - \theta)g(Z_1)]$  on
square-integrable functions $g$ (that is,  $\Exp[g(Z_1)^2]<\infty$).

\begin{thm} \label{theorem: degenerate asymptotics} If the kernel $k$
  of the statistic $U_n$ from~(\ref{eq:Un}) has variance
  $\sigma_m^2 < \infty$ and $\sigma_1^2 = 0$, then
\begin{align*}
	n(U_n - \theta) \overset{d}{\to} {m \choose 2}\sum_{i=1}^\infty \lambda_i (\chi_{1i}^2 - 1)
\end{align*}
where $\chi_{11}^2,\chi_{12}^2,\dots$ are i.i.d.~$\chi^2_1$ random
variables, and the $\lambda_i$'s are the eigenvalues, taken with
multiplicity, associated to a system of orthonormal eigenfunctions of
  the operator $A_k$.
\end{thm}

We will use Theorem \ref{theorem: degenerate asymptotics} in Section
\ref{section: asymptotics under the null} to find the asymptotic
distribution of $t^*$ under the null hypothesis of independence.

\section{Degeneracy of the sign covariance}
\label{sec:degeneracy}

Let $Z_i=(X_i,Y_i)$ for $i = 1,2,\dots,n$ be an i.i.d.\ sequence
comprising copies of a random vector $(X,Y)$ with values in $\mathbb{R}^2$.
Let $t^*$ be the (empirical) Bergsma-Dassios sign covariance for this
sample.  We begin our study of the asymptotic properties of the
U-statistic $t^*$ by studying its degeneracy.  Our first observation
is that $t^*$ is degenerate when $X$ and $Y$ are independent, denoted
$X\perp\!\!\!\perp Y$.  Next, in a particular setting that has $(X,Y)$
continuously distributed, we are able to show that $t^*$ is degenerate
only if $X\perp\!\!\!\perp Y$.

The statistic $t^*$ has the kernel $h$ from \eqref{eq: kernel}, which
has $m=4$ arguments.  Specializing the definitions from~(\ref{eq:k_i})
and~(\ref{eq:sigma_i^2}) to the present setting, we may define
functions $h_1,...,h_4$ with variances $\sigma_1^2,...,\sigma_4^2$.
The kernel $h$ is a bounded function and thus $\sigma_4^2<\infty$.
Hence, Theorem~\ref{theorem: normal asymptotics} applies and yields
the following result.

\begin{cor}\label{cor: asymptotics under dependence}
  As $n\to\infty$, the sign covariance converges to a normal limit,
  namely,
  \begin{align*}
    \sqrt{n}(t^* - \tau^*) \overset{d}{\to} N(0, 16\sigma_1^2).
  \end{align*}
\end{cor}

The result just stated provides a non-trivial distributional
approximation to $t^*$ only if $\sigma_1^2>0$.  The following lemma
observes that this fails to be the case under the null hypothesis of
independence, under which $t^*$ is a degenerate U-statistic.  The
proof of the lemma as well as the proofs of all other results in this
section are deferred to Appendix~\ref{section: proofs-degen}.

\begin{lem}\label{lemma: degeneracy of h_1}
  If $X\perp\!\!\!\perp Y$ then $\sigma_1^2 = \Var[h_1(X_1,Y_1)] =0$ so that
  $h_1(X_1,Y_1)$ is a degenerate random variable.
\end{lem}

According to Lemma~\ref{lemma: degeneracy of h_1} and Theorem
\ref{theorem: normal asymptotics}, if $X\perp\!\!\!\perp Y$ we have
$\sqrt{n}t^*\overset{p}{\to} 0$ because $\Exp[t^*] =\tau^*= 0$ under
independence \citep{BergsmaDassios14}.  We thus need to appeal to
Theorem \ref{theorem: degenerate asymptotics} to find a non-degenerate
asymptotic distribution for $t^*$ when $X\perp\!\!\!\perp Y$.  This is the topic
of Section~\ref{section: asymptotics under the null}.

\begin{rem} \label{remark: sigma values} In the continuous case with
  $X\perp\!\!\!\perp Y$, it is possible to compute all of the variances
  $\sigma_1^2,...,\sigma_4^2$ exactly.  We report these values to be
  \begin{align*}
    \sigma_1^2 = 0,\ \ \ \sigma_2^2 = \frac{1}{225}, \ \ \ \sigma_3^2 = \frac{8}{225}, \ \ \ \sigma_4^2 = \frac{50}{225}.
  \end{align*}
  The fact that $\sigma_1^2=0$ was shown in generality in
  Lemma~\ref{lemma: degeneracy of h_1}.  The value of $\sigma_2^2$ can
  be computed as the sum of the squared eigenvalues of $h_2$ which are
  derived in the proof of Theorem \ref{theorem: continuous asymptotic
    distribution}; in particular, we have that
  $\sigma_2^2 = \sum_{i=1}^\infty\sum_{j=1}^\infty \frac{6^2}{\pi^8
    i^4j^4} = \frac{1}{225}$. Finally, $\sigma_3^2$ and $\sigma_4^2$
  can be computed explicitly from the representation of $h$ in Lemma
  \ref{lemma: kernel values}, this computation is trivial for
  $\sigma_4^2$ but quite lengthy for $\sigma_3^2$ and thus is omitted.
\end{rem}

Next, we turn our attention to the case that $X\not\perp\!\!\!\perp Y$ and
$(X,Y)$ are generated from a continuous distribution on
$\mathbb{R}^2$.  In this case we find $t^*$ to be non-degenerate.

\begin{thm}\label{theorem: nondegeneracy of h_1}
  Suppose $(X,Y)$ has a bivariate continuous distribution with a
  continuous density function $f$ with support
  $\ol{f^{-1}((0, \infty))} = [a,b]\times [c,d]$, where
  $-\infty\le a<b\le \infty$ and $-\infty\le c<d\le \infty$.  If $X$
  and $Y$ are dependent, then $\sigma_1^2 = \Var[h_1(X_1,Y_1)] > 0$.
\end{thm}

In the setting of Theorem~\ref{theorem: nondegeneracy of h_1}, we thus
have that $t^* = N(\tau^*, 16 \sigma_1^2/n) + o_p(n^{-1/2})$.

\begin{exmp} \label{example: variance of h1} To gain intuition for the
  magnitude of the asymptotic variance $16\sigma_1^2$, we use Monte
  Carlo integration to compute $16\sigma_1^2$ in the case that $(X,Y)$
  follow a bivariate normal distribution.  Since
  $\sigma_1^2$ is an even function of the correlation $\rho$ of a
  bivariate normal distribution, we consider $\rho \in [0,1]$.  In
  particular, we perform this computation letting $\rho$ take on 20,
  evenly spaced, values between $0$ and $1$.  The results of this
  computation are shown in Figure \ref{figure: variance of h1}. The
  figure shows that the asymptotic variance $16 \sigma_1^2$ gradually
  increases with larger values of $\rho$, peaking at $16 \sigma_1^2
  \approx 0.14$ when $\rho \approx 0.74$, and then decreases to 0 as
  the correlation further approaches 1. Note that a value of
  $\sigma_1^2 = 0$ when $\rho=1$ does not contradict
  Theorem~\ref{theorem: nondegeneracy of h_1} as, in this case, the
  joint distribution of $(X,Y)$ is not continuous. The shape of the curve in Figure
  \ref{figure: variance of h1} can be partially explained by the fact that
  $\sigma_1^2 \leq \sigma_4^2 = \Var(h(Z_1,...,Z_4)) = (\tau^* + 1/3)(1 - (\tau^* + 1/3))$.
  For instance, note that $(\tau^* + 1/3)(1 - (\tau^* + 1/3))$ equals 0 when $\tau^* = 2/3$
  (in which case the correlation can be seen to be 1 or -1), and is maximized at
  $\tau^* = 1/6 \approx .167$ which corresponds a correlation of approximately 0.7 (see
  Figure \ref{fig:normal-taustar}).
\end{exmp}

\begin{figure}
        \centering
        \includegraphics[width=.75\textwidth]{h1_variance.pdf}
        \caption{Monte Carlo approximations to the values of
          $16\sigma_1^2$ for bivariate normal distribution with
          different correlations.}
        \label{figure: variance of h1}
\end{figure}

\section{Asymptotics under the null hypothesis of
  independence} \label{section: asymptotics under the null}

As in the previous section, let $t^*$ be the empirical sign covariance
for an i.i.d.~sample $(X_i,Y_i)$, $i = 1,2,\dots,n$, with values in
$\mathbb{R}^2$.  Throughout this section, we assume the $(X_i,Y_i)$ to be
independent copies of a random vector $(X,Y)$ with $X\perp\!\!\!\perp Y$, so
that $t^*$ is degenerate (Lemma~\ref{lemma: degeneracy of h_1}).  We
thus need to appeal to Theorem \ref{theorem: degenerate asymptotics}
to find a non-degenerate asymptotic distribution for $t^*$.  Since
$\Exp[t^*] =\tau^*= 0$ under independence, we are led to the problem
of determining the eigenvalues of the operator
$A_h:g(\cdot) \mapsto \Exp[h_2(\cdot, Z_1)g(Z_1)]$.

A key observation is that under independence $A_h$ is a tensor product
of operators because the function $h_2$ admits the following
factorization, which along with all other results in this section is
proved in Appendix~\ref{section: proofs}.

\begin{lem}\label{lemma: product formula for h_2}
  If $X\perp\!\!\!\perp Y$ then
  \[
    h_2((x_1,y_1),(x_2,y_2)) = \frac{2}{3}\,g_{X}(x_1,x_2)
    \,g_{Y}(y_1,y_2),
  \]
  where $g_{X}(x_1,x_2) = \Exp[a(x_1,x_2,X_3,X_4)]$
  and $g_{Y}(y_1,y_2) = \Exp[a(y_1,y_2,Y_3,Y_4)]$.
\end{lem}

The function $g_X$ (and similarly $g_Y$) takes the form
\begin{align}
  \label{eq:prob-form-of-g}
	 g_{X}(x_1,x_2)
 	&\ = P(x_1, X_3 < x_2, X_4) + P(x_1, X_3 > x_2, X_4) \\
  \nonumber
	&\hspace{6mm} - P( x_1, x_2 < X_3, X_4) - P(x_1, x_2 > X_3, X_4).
\end{align}

By Lemma~\ref{lemma: product formula for h_2}, $A_h=A_{g_X}\otimes
A_{g_Y}$ and thus the spectrum of $A_h$ is the product of the spectra
of $A_{g_X}$ and $A_{g_Y}$.  We record the general version of this
fact in the next lemma.  Here, eigenvalues are always repeated
according to their multiplicity, and we let
$\mathbb{N}_+=\{1,2,\dots\}$.

\begin{lem} \label{lemma: h_values} Let $g_1$ and $g_2$ be symmetric
  real-valued functions with $\Exp[g_1(X_1,X_2)] = \Exp [g_2(Y_1,Y_2)]
  = 0$ and $\Exp[g_1(X_1,X_2)^2], \Exp [g_2(Y_1,Y_2)^2] <\infty$.  For
  $i=1,2$, let $\lambda_{i,j}$, $j\in\mathbb{N}_+$, be the nonzero eigenvalues
  of $A_{g_i}$.

  Then the products $\lambda_{1,j_1}\lambda_{2,j_2}$, $j_1,j_2\in\mathbb{N}_+$,
  are the nonzero eigenvalues of $A_k$ for $k((x_1,y_1),(x_2,y_2)) :=
  g_1(x_1,x_2)g_2(y_1,y_2)$.
\end{lem}

In the sequel, we use the factorization results from
Lemmas~\ref{lemma: product formula for h_2} and~\ref{lemma: h_values}
to obtain the asymptotic distribution of $t^*$ when $X$ and $Y$ are
continuous (Section~\ref{subsection: continuos}), and when $X$ and $Y$
are discrete with finite support (Section~\ref{subsection: discrete}).
A straightforward extension covers the mixed continuous and discrete
case (Section~\ref{subsection: discrete}).

\subsection{Continuous variables}\label{subsection: continuos}

Suppose now that $X\perp\!\!\!\perp Y$ with $X$ and $Y$ following continuous
marginal distributions.  Since $h((X_1,Y_1),...,(X_4,Y_4))$ depends
only on the joint ranks of $(X_1,Y_1),...,(X_4,Y_4)$, it follows that
$\tau^*$ (and $t^*$) are invariant to monotonically increasing
transformations of the marginals of $(X,Y)$.  As such we may, and
will, assume that $X$ and $Y$ are i.i.d.~$\text{Uniform}(0,1)$.  Then
$(X,Y)$ is uniform on the unit square $(0,1)\times (0,1)$.  In this
case the factorization described in Lemma \ref{lemma: product formula
  for h_2} has a particularly nice form.

\begin{lem} \label{lemma: continuous factorization}
  If $X,Y\overset{i.i.d.}{\sim}\text{Uniform}(0,1)$, then for $(x_1,y_1),(x_2,y_2) \in (0,1)^2$,
  \begin{align*}
    h_2((x_1, y_1), (x_2, y_2)) &= 6\, c(x_1,x_2)\,c(y_1,y_2)
\end{align*}
where
\[
  c(x_1,x_2)=\frac{1}{2} x_1^2+ \frac{1}{2}x_2^2- x_1\vee
  x_2+\frac{1}{3}
\]
and $x_1\vee x_2:=\max\{x_1,x_2\}$.
\end{lem}

Somewhat surprisingly, the function $c$ corresponds to the kernel of
the well studied Cram\'{e}r-von Mises statistic.  Leveraging the fact
that the eigenvalues of $A_c$ are already known, we are now able to
derive the asymptotic distribution of $t^*$.

\begin{thm} \label{theorem: continuous asymptotic distribution} If $X$
  and $Y$ are independent continuous random variables, then
\begin{align*}
	nt^* &\overset{d}{\to} \frac{36}{\pi^4}\sum_{i=1}^\infty\sum_{j=1}^\infty \frac{1}{i^2j^2}(\chi_{1,ij}^2-1)

\end{align*}
where $\{\chi^2_{1,ij} : i,j \in \mathbb{N}_+ \}$ is a collection of i.i.d.\
$\chi^2_1$ random variables.
\end{thm}

Remarkably the asymptotic distribution just given is simply a scale
multiple of the asymptotic distribution of the U-statistic for
Hoeffding's $D$ where
\[
  D=\iint (F_{X,Y}(x,y)-F_X(x)F_Y(y))^2 \ \text{d}F_{X,Y}(x,y);
\]
see \cite{Hoeffding48}.  When $(X,Y)$ has a continuous joint
distribution, it is readily seen that $D=0$ if and only if
$X\perp\!\!\!\perp Y$.  However, this may fail in non-continuous cases.

\subsection{Discrete variables} \label{subsection: discrete}

We now treat the case where $X$ and $Y$ are independent discrete
random variables with finite supports.  Unlike in the continuous case,
the asymptotic distribution of $t^*$ then depends on how $X$ and $Y$
distribute their probability mass marginally. In practical applications these
marginal probabilities must be estimated before using our limit theorem.

In order to present the result, we associate a matrix to a discrete
random variable as follows.  Let $U$ be a random variable with finite
support $\{u_1,\dots,u_r\}$, cumulative distribution function $F_U$
and probability mass function $p_U$.  We then define $R^U$
to be the $r\times r$ symmetric matrix whose $(i,j)$-th entry is
\begin{align} \label{equation: R matrix form}
	R^U_{ij} =&
	\ \sqrt{p_U(u_i)p_U(u_j)} \Bigg\{ \bigg[(F_U(u_i\wedge u_j) - p_U(u_i\wedge u_j))^2 + (1-F_U(u_i\vee u_j))^2\bigg] \\
	&- I(u_i\not=u_j)\bigg[F_U(u_i\wedge u_j)(1-F_U(u_i\wedge u_j))
	+ \sum_{u_i\wedge u_j < u_\ell  < u_i\vee u_j} p_U(u_\ell)(1-F_U(u_\ell))\bigg] \nonumber
	 \Bigg\}.
\end{align}

\begin{thm}\label{theorem: discrete asymptotic distribution}
  Let $X$ and $Y$ be independent discrete random variables with
  finite supports of size $r$ and $s$, respectively.  Let
  $\lambda_1^X,\ldots,\lambda_r^X$ be the eigenvalues of $R^X$, and
  let $\lambda^Y_1,\ldots,\lambda^Y_s$ be the eigenvalues of $R^Y$.
  Then
 \begin{align*}
	nt^* &\overset{d}{\to} 4\sum_{i=1}^r\sum_{j=1}^s \lambda_i^X \lambda_j^Y (\chi_{1,ij}^2-1)
\end{align*}
where $\{\chi^2_{1,ij} : i \leq r,~j\leq s \}$ is a collection of $rs$ i.i.d.\ $\chi_1^2$ random variables.
\end{thm}

In the special case that $X$ and $Y$ are Bernoulli random variables,
the asymptotic distribution can be presented in simple form.

\begin{exmp}
  If $X\sim\text{Bernoulli}(p)$ for $p\in(0,1)$, then
  \begin{align*}
    R^X = \begin{pmatrix}
      p^2(1-p) & -(p(1-p))^{3/2} \\
      -(p(1-p))^{3/2} & p(1-p)^2
    \end{pmatrix}
  \end{align*}
  has rank one and its nonzero eigenvalue is $p(1-p)$.  It follows
  that if $Y$ is a second independent random variable with
  $Y\sim\text{Bernoulli}(q)$ for $q\in(0,1)$, then
  \[
    nt^* \overset{d}{\to} 4pq(1-p)(1-q)(\chi^2_1 - 1).
  \]
  So, $t^*$ can be centered and scaled to become asymptotically
  chi-square.
\end{exmp}

\begin{exmp}
  For a ternary random variable $X$ with $P(X=1)=p_1$, $P(X=2)=p_2$
  and $P(X=3)=p_3=1-p_1-p_2$, we have
  \begin{align*}
    R^X =  \begin{pmatrix}
      p_1(1-p_1)^2
      & -\sqrt{p_1p_2}\left[ p_1(1-p_1)-p_3^2\right]
      & -\sqrt{p_1p_3}\left[ p_3(1-p_3)+p_1p_2\right]\\
      .
      & p_2\left( p_1^2+p_3^2\right)
      & -\sqrt{p_2p_3}\left[p_3 (1 -p_3)- p_1^2 \right]\\
      .
      & .
      &p_3(1-p_3)^2\\
    \end{pmatrix},
  \end{align*}
  where we show only the upper half of the symmetric matrix.  No
  simple formula seems to be available to determine the eigenvalues of
  $R^X$ in this case, but the eigenvalues can readily be computed
  numerically for any (possibly estimated) values of $p_1$ and $p_2$.
\end{exmp}

Finally,  if $X$ is discrete with finite support and $Y$ is
continuous, then a simple extension of Theorems \ref{theorem:
  continuous asymptotic distribution} and \ref{theorem: discrete
  asymptotic distribution} gives the following result.

\begin{cor}\label{theorem: mixed asymptotic distribution}
  Let $X$ and $Y$ be independent random variables, where $X$ has
  finite support of size $r$ and $Y$ is continuous. Let $\lambda_1,...,\lambda_r$
  be the eigenvalues of $R^X$. Then
 \begin{align*}
   nt^* &\overset{d}{\to} \frac{12}{\pi^2}\sum_{i=1}^r\sum_{j=1}^\infty \frac{\lambda_i}{j^2} (\chi_{1,ij}^2-1)
\end{align*}
where $\{\chi^2_{1,ij} : i \leq r,~j \in\mathbb{N}_+ \}$ is a collection of
i.i.d.\ $\chi^2_1$ random variables.
\end{cor}

\section{Simulations}\label{section: simulations}

The results from Section~\ref{section: asymptotics under the null} can
be used to form asymptotic tests of independence, and we now explore
which sample sizes are needed for the asymptotic approximations to be
accurate.  As a test based on $t^*$ has asymptotic power against all
alternatives to independence, it is also of interest to make
comparisons against other tests known to be (most) powerful for
particular settings and alternatives.  Finally, we demonstrate how the
results of Section~\ref{sec:degeneracy} can be used for sample size
computations.  Code for performing asymptotic tests of independence
has been
incorporated in the TauStar\footnote{See
 \url{https://cran.r-project.org/web/packages/TauStar/index.html}} R
package available on CRAN, the Comprehensive R Archive Network
\citep{R, taustar-package}.

\subsection{Empirical convergence to the asymptotic
  distribution} \label{subsection: empirical convergence}

Let $t^*$ be computed from a sample of size $n$ drawn from the joint
distribution of a bivariate random vector $(X,Y)$ with $X\perp\!\!\!\perp Y$.
Since $t^*$ only depends on ranks, its distribution does not change
when applying monotonically increasing marginal transformations to $X$
and $Y$.  When $X$ and $Y$ both have continuous distributions, we may
thus transform their distributions to $N(0,1)$ without changing the
distribution of $t^*$.  When one or both of $X$ and $Y$ are discrete
however, the distribution of $t^*$ depends on how $X$ and $Y$
distribute their probability mass making it impossible to provide an
exhaustive empirical study of convergence properties.  Instead we will
consider selected examples.  Specifically, we consider the following cases:
\begin{enumerate}[(i)]
	\item The continuous case with $X,Y\sim N(0,1)$.
	\item A discrete case with $P(X=i) = 1/10$ for $1\leq i\leq
          10$, and $P(Y=i) \propto 2^{-i}$ for $1\leq i\leq 12$.
	\item A mixed case with $X\sim N(0,1)$ and $P(Y = i) = 1/5$
          for $1\leq i\leq 5$.
\end{enumerate}
In each setting we compute, for different sample sizes $n$, a kernel
density estimate for the distribution of $t^*$ and plot it alongside
the asymptotic density.  The resulting plots are shown in Figure
\ref{figure: density convergence}, which demonstrates that the
asymptotic and finite-sample distributions are in close agreement
already when $n=80$.  While we present only one example each for the
discrete and mixed cases we found similar results when simulating with
many other choices of distributions.

\begin{figure}
        \centering
         \subfloat[Continuous case\label{figure: continuous density}]{%
		\includegraphics[width=.33\textwidth]{cont_dens_conv}
	}
	\subfloat[Discrete case\label{figure: discrete density}]{%
		\includegraphics[width=.33\textwidth]{dis_dens_conv}
	}
	\subfloat[Mixed case\label{figure: mixed density}]{%
		\includegraphics[width=.33\textwidth]{mixed_dens_conv}
	}
        \caption{Kernel density estimates from 10,000 simulated values
          of $nt^*$ at sample sizes
          $n \in \{10, 15, 20, 25, 30, 40, 50, 60, 70, 80\}$; smaller
          values of $n$ are shown in lighter color.  The plots also
          show the density of the asymptotic distributions from
          Section \ref{section: asymptotics under the null} in dashed
          blue line.}
        \label{figure: density convergence}
\end{figure}

\begin{rem}
  Computing the asymptotic densities shown in Figure \ref{figure:
    density convergence} is non-trivial and requires the numerical
  inversion of the characteristic function for the asymptotic
  distributions.  To perform this numerical inversion we use the
  techniques described in Section 7 of \citet*{BlumEtAl61}; these
  computations are done automatically in the aforementioned TauStar
  package for R.
\end{rem}

\subsection{Power comparisons} \label{subsection: power comparisons}

We explore the power of an asymptotic test based on $t^*$ in six
cases:
\begin{enumerate}[(i)]
\item First, we take $(X,Y)$ as bivariate normal with correlation
  $\rho \in\{0, .1, .2, \dots, 1\}$; the distribution of $t^*$ then does
  not depend on the means and variances which may thus be set to zero
  and one, respectively.  We compare the test based on $t^*$ to the
  two-sided test based on the standard Pearson
  correlation $\hat{\rho}$.  We implement the latter test using the
  fact that $\hat{\rho}\sqrt{(n-2)/(1-\hat{\rho}^2)}$ has a
  $t$-distribution with $n-2$ degrees of freedom.

\item Next, we consider three discrete cases all of which have $(X,Y)$
  taking values in the grid $\{1,2,\dots,5\}^2$.  In each of these cases
  we compare our test to the chi-square test of independence.
  \begin{enumerate}[(a)]
  \item In the first discrete case, $(X,Y)$ follows a mixture between
    the uniform distribution on $\{(1,1),(2,2),\dots,(5,5)\}$ and the
    uniform distribution on $\{1,\dots,5\}^2$, illustrated in Figures
    \ref{figure: dg ordinal} and \ref{figure: dg uniform},
    respectively.  We let the mixture weight $p$ for the former
    component range through the set $\{0,.1,\dots,1\}$.
  \item The second discrete case is analogous but a mixture between
    the distributions from Figures \ref{figure: dg nonordinal} and
    \ref{figure: dg uniform}.
  \item The third discrete case is as the previous two but mixes the
    distributions from Figures \ref{figure: dg L shape} and
    \ref{figure: dg uniform}.
  \end{enumerate}

\begin{figure}
        \centering
         \subfloat[\label{figure: dg ordinal}]{%
		\begin{tikzpicture}
			\fill[black!20!white] (0,2) rectangle (.5,2.5);
			\fill[black!20!white] (.5,1.5) rectangle (1,2);
			\fill[black!20!white] (1,1) rectangle (1.5,1.5);
			\fill[black!20!white] (1.5,.5) rectangle (2,1);
			\fill[black!20!white] (2,0) rectangle (2.5,.5);
			\draw[step=0.5cm,black,very thin] (0,0) grid (2.5,2.5);
		\end{tikzpicture}
	}
	\hspace{10mm}
	\subfloat[\label{figure: dg nonordinal}]{%
		\begin{tikzpicture}
			\fill[black!20!white] (0,2) rectangle (.5,2.5);
			\fill[black!20!white] (2,1.5) rectangle (2.5,2);
			\fill[black!20!white] (1,1) rectangle (1.5,1.5);
			\fill[black!20!white] (1.5,.5) rectangle (2,1);
			\fill[black!20!white] (.5,0) rectangle (1,.5);
			\draw[step=0.5cm,black,very thin] (0,0) grid (2.5,2.5);
		\end{tikzpicture}
	}
	\hspace{10mm}
        	\subfloat[\label{figure: dg L shape}]{%
		\begin{tikzpicture}
			\fill[black!20!white] (0,0) rectangle (2.5,.5);
			\fill[black!20!white] (0,0) rectangle (.5,2.5);
			\draw[step=0.5cm,black,very thin] (0,0) grid (2.5,2.5);
		\end{tikzpicture}
	}
	\hspace{10mm}
	\subfloat[\label{figure: dg uniform}]{%
		\begin{tikzpicture}
			\fill[black!20!white] (0,0) rectangle (2.5,2.5);
			\draw[step=0.5cm,black,very thin] (0,0) grid (2.5,2.5);
		\end{tikzpicture}
	}
        \caption{Visualization of where probability mass is placed for
          different discrete distributions on $\{1,...,5\}^2$. In each
          case the distribution is uniform over the gray squares, and
          zero probability is assigned to the white squares.} \label{figure: discrete grids}
\end{figure}

\item Finally, we experiment with two mixed cases in which the
  distribution of $X$ is
  discrete and the conditional distribution of $Y$ given $X=x$ is
  a normal distribution $N(\mu_x,1)$.
  \begin{enumerate}[(a)]
  \item The first mixed case has $X\sim\text{Bernoulli}(.3)$,
    $\mu_x=0$ when $x=1$, and $\mu_x=\mu$ when $x=0$.  Here, we let
    the mean difference $\mu$ range through the set $\{0, 1/6,
    2/6,..., 9/6\}$, and each setting we compare against the
    two-sample $t$-test.
  \item In the second case $X\sim\text{Uniform}(\{1,...,6\})$, and $Y$
    has conditional mean $\mu_x$ is zero when $x$ is odd and equal to
    $\mu \in \{0, 1/6, 2/6,..., 9/6\}$ when $x$ is even.  Here, we
    compare against a bootstrapped permutation test using the distance
    covariance statistic of \cite{SzekelyEtAl07}, known to be
    consistent for independence, using the Energy R package
    \cite{RizzoSzekely14}.
  \end{enumerate}
\end{enumerate}

The simulation results are presented in Figure \ref{figure: power
  simulations}.  Surprisingly, the $t^*$ test has competitive power in
cases (i) and (iii)(a) where the alternative tests are known to be
most powerful given the distributional assumption of normality. For
the jointly discrete cases, we observe that the chi-square test of
independence has essentially equal power in case (ii)(a),
significantly higher power in case (ii)(b), and significantly lower
power in case (ii)(c). The lack of power in case (ii)(b) is not
surprising as the $t^*$ statistic is ordinal in nature and the
dependence in the distribution from case (ii)(b) was designed to be
non-ordinal.  The ordinal nature of $t^*$ also explains the
significant gains in case (ii)(c).  Hence, it would seem that the
$t^*$ test for jointly discrete data can offer substantial
improvements in power over the chi-square test if an ordinal
dependence relationship is suspected in the data. Finally, case
(iii)(b) suggests that there are cases in which $t^*$ may provide
higher power than the distance covariance.

\begin{figure}
        \centering
         \subfloat[Case (i)\label{figure: case (i)}]{%
         		\includegraphics[width=.35\textwidth]{cont_levels_plot}
	}
	\subfloat[Case (ii)(a)\label{figure: case (ii)(a)}]{%
         		\includegraphics[width=.35\textwidth]{dis_ord_levels_plot}
	}
	~ \\
        	\subfloat[Case (ii)(b)\label{figure: case (ii)(b)}]{%
         		\includegraphics[width=.35\textwidth]{dis_nonord_levels_plot}
	}
	\subfloat[Case (ii)(c)\label{figure: case (ii)(c)}]{%
         		\includegraphics[width=.35\textwidth]{dis_L_shape_levels_plot}
	}
	~ \\
	\subfloat[Case (iii)(a)\label{figure: case (iii)(a)}]{%
         		\includegraphics[width=.35\textwidth]{mixed_binary_levels_plot}
	}
	\subfloat[Case (iii)(b)\label{figure: case (iii)(b)}]{%
         		\includegraphics[width=.35\textwidth]{mixed_zig_levels_plot}
	}
        \caption{Simulated power of the $t^*$ asymptotic test (in
          solid line) and the power of the competing test (in dashed
          line). In each case we use a level of 0.05, which is
          displayed as a horizontal dashed black line.} \label{figure: power
          simulations}
\end{figure}

\subsection{Sample size calculations} \label{subsection: computing power}

Focusing on the continuous case, consider an asymptotic level $\alpha$
test of the null hypothesis of $\tau^*=0$ (i.e., independence) that
compares the statistic $t^*$ to a critical value $c_{\alpha}$ derived
from the asymptotic distribution from Theorem \ref{theorem: continuous
  asymptotic distribution}.  Suppose we would like to determine
the minimum sample size $n_\beta$ needed for a power of at least
$\beta$ under an alternative that has the two considered variables $X$
and $Y$ dependent, so that $\tau^*>0$.  If the sample is drawn from a
joint distribution for $(X,Y)$ that satisfies the conditions of
Corollary \ref{cor: asymptotics under dependence}, and if $\sigma_1^2$
is known to be no larger than the quantity $\bar\sigma_1^2$, then
Corollary \ref{cor: asymptotics under dependence} implies that for any
$x \leq \tau^*$,
\begin{align*}
  P(t^*\leq x) &= P\left(\sqrt{n}(t^*-\tau^*) \leq \sqrt{n}(x-\tau^*)\right) \\
               &\approx P\left(N(0,16\, \sigma_1^2) \leq
                 \sqrt{n}(x-\tau^*)\right) \\
               &\leq P\left(N\left(\tau^*, {16\, \bar\sigma_1^2
                 }/{n}\right) \leq x\right).
\end{align*}
This result can be used to find an asymptotically valid upper bound
$\bar n_\beta$ on $n_\beta$, by letting $\bar n_\beta$ be the smallest
positive integer such that
$c_{\alpha}/\bar n_\beta \leq \tau^*$ and
$P(N(\tau^*, {16\, \bar\sigma_1^2}/{\bar n_\beta}) \leq
c_\alpha/\bar n_\beta) \leq 1-\beta$. Finding this number
$\bar n_\beta$ can be accomplished in an iterative fashion.

The remaining difficulty in such an asymptotic sample size calculation
is finding a suitable upper bound $\bar\sigma_1^2$ for the unknown
variance $\sigma_1^2$. A crude but universally valid upper bound for
$\sigma_1^2$ can be obtained from Lemma \ref{lemma: kernel values},
which implies that $h_1(X,Y)$ takes values in the interval
$[-1/3,2/3]$ and thus $\sigma_1^2 \leq 1/4$. When $(X,Y)$ is bivariate normal,
an approximately valid upper bound of $\sigma_1^2$ is given by
$\sigma_1^2 \leq 0.14/16 = 0.00875$ (see Example \ref{example: variance
  of h1}). Figure \ref{figure: power upper bound} plots the
upper bound for the minimum sample size needed to achieve various powers when
bounding $\sigma_1^2$ by $1/4$ and $0.00875$ respectively and sampling
from a bivariate normal distribution with correlation 0.6. From the figure, we see that the $1/4$
bound leads to very conservative sample sizes while the $0.00875$ bound
results in values that much more closely adhere to the empirical
truth.  In general, overestimation of $\sigma_1^2$ is advisable as
small values of $\bar\sigma_1^2$ may lead to consideration of sample
sizes that are too small for asymptotic approximations to be reflective
of the actual finite-sample behavior of the test.

\begin{figure}
        \centering
        \includegraphics[width=.75\textwidth]{power_upper_bound.pdf}
        \caption{Minimum sample size ($n_{\beta}$) needed to a achieve a desired
          power $\beta$ at level $0.05$.  Simulations for bivariate normal
          data with correlation $0.6$ were used to compute an estimate
          of $n_{\beta}$ (blue line with dots).  These are
          compared to two asymptotic upper bounds for $n_{\beta}$
          using the bound $\sigma_1^2\leq 1/4$ (red line) and the
          bound $\sigma_1^2 \leq 0.00875$ (dashed orange line).  The
          sample size is presented with a log scaling.} \label{figure:
          power upper bound}
\end{figure}

\section{Discussion}\label{section: discussion}

The sign covariance $\tau^*$ of \cite{BergsmaDassios14} has the
intriguing property of being zero if and only if the considered pair
of random variables is independent, assuming that the random variables
follow a distribution that is continuous, discrete or a mixture of
such distributions.  Under these mild conditions, testing the
hypothesis that $\tau^*=0$ thus allows one to consistently assess
(in-)dependence.  With the aim of simplifying the implementation of
such independence tests, we have given a comprehensive study of the
asymptotic properties of $t^*$, the natural U-statistic for $\tau^*$.
The asymptotic distribution of $t^*$, especially as described in
Section \ref{subsection: continuos}, is seen to be connected in
interesting ways to the asymptotic distribution of Hoeffding's $D$,
and the Cram\'{e}r-von Mises statistic.

One limitation of our work is that we did not consider asymptotic
distributions under local alternatives to independence.  The reason is
that these would be distributions of weighted sums of non-central
chi-square random variables, which seem difficult to use in numerical
computations for assessment of power or sample size calculation.

While we have a complete understanding of the asymptotics of $t^*$
under fairly weak distributional assumptions---we covered continuous
and discrete cases, it remains to be seen if the large-sample
distribution of $t^*$ can be obtained without any such assumptions.
However, as noted above, it is also not yet known if the property that
$\tau^*=0$ only under independence holds for distributions that
are not continuous, discrete or a mixture of two such distributions.


\section{Introduction}
Throughout this paper $G = (V, E)$ is a  connected finite simple undirected graph with vertex set $V = V(G)$ and edge set $E =E(G)$.  Let $|G|$ or $|V|$ denote the cardinality of $V$.
For $S \subseteq V(G)$ and  $F \subseteq E(G)$,   we use $G[S]$ for the subgraph of $G$ induced by $S$, $G - S$ for the subgraph induced by $V(G) - S$ and $G - F$ for the subgraph of G obtained by deleting $F$.  Let $w(G - S)$ be the number of components of $G - S$, and  $S$ be a cut set if $w(G - S) \geq 2$.  For a vertex $v \in V(G)$, the neighborhood of  $v $ is the set $N(v) = N_G(v) = \{w \in V(G), vw \in E(G)\}$, and $d_G(v)$ (or briefly $d(v)$)  denotes the degree of $v$ with $ d_{G}(v) = |N(v)|$.
 We use $G \cong H$ to denote that $G$ is isomorphic to $H$ and $G \not\cong H$ to denote  that $G$ is not isomorphic to $H$. Let $K_n, P_n$ and $S_n$ denote the clique, the path and the star on $n$ vertices, respectively. In particular, we say $K_n$ is a $k$-clique for $n = k$.

In the 1980s, Narumi and Katayama [7] considered the product
$$\begin{array}{rcl}
NK = \prod_{v \in V(G)} d(v)
\end{array}$$
which is  the "Narumi-Katayama index". And recently, Todeschini and Gutman et al [4, 10, 11]  studied the first and second multiplicative Zagreb indices defined as follow:
$$\begin{array}{rcl}
\prod_1(G) &=& \prod_{v \in V(G)}d(v)^{2},\\
\prod_2(G) &=& \prod_{uv \in E(G)} d(u)d(v).
\end{array}$$
Obviously, the first multiplicative Zagreb index is just the square of the NK index. Gutman [4] in 2011 characterized the multipilicative Zagreb indices for trees and determined the unique trees that obtained maximum and minimun values for $\prod_1(G)$ and $\prod_2(G)$, respectively.

\vskip 2mm {\bf Theorem 1 (Gutman 2011) } \emph{ Let $n \geq 5$ and  $T_n$ be any tree with $n$ vertices, then
$$\begin{array}{rcl}
(i) \prod_1(S_n) \leq \prod_1(T_n) \leq \prod_1(P_n);\\
 (ii)\prod_2(P_n) \leq \prod_2(T_n) \leq \prod_2(S_n).
\end{array}$$
  }

In this paper, we consider the first generalized multiplicative Zagreb index defined in (1) below and the second multiplicative Zagreb index:
for any real number $c>0$,
$$\begin{array}{rcl}
(1)\prod_{1,c}(G) &=& \prod_{v \in V(G)}d(v)^{c};\\
(2)\prod_2(G) &=& \prod_{uv \in E(G)} d(u)d(v).
\end{array}$$
Eventually, for $c = 1, 2$,  (1) is just the NK index and the first multuplicative Zagreb, respectively. For (2), it is easy to see that $\prod_2(G) = \prod_{v\in V(G)}d(v)^{d(v)}$. Also we will find the bounds of the values of   $\prod_{1,c}(G)$, $\prod_{2}(G)$ for $k$-trees, respectively, and determine the extremal graphs which attain the bounds. Our main results are as follows:

{\bf Theorem 2 }  \emph
Let $T_n^k$ be a $k$-tree on $n \geq k$ vertices, then
 $$\begin{array}{rcl}
  \prod_{1,c}(S_{k,n-k}) \leq \prod_{1,c} (T_n^k) \leq \prod_{1,c}(P_n^k),\end{array}$$
the left-side  and the right-side equalities are reached if and only if $T_n^k \cong S_{k,n-k}$ and $T_n^k \cong P_n^k$, respectively.

\vskip 2mm
{\bf Theorem 3 }  \emph
Let $T_n^k$ be a $k$-tree on $n \geq k$ vertices, then
 $$\begin{array}{rcl}
  \prod_2(P_n^k) \leq \prod_2 (T_n^k) \leq \prod_2(S_{k,n-k}),\end{array}$$
 the left-side  and the right-side equalities are reached if and only if $T_n^k \cong P_n^k$ and $T_n^k \cong S_{k,n-k}$, respectively.

\section{Preliminary}
It is commonly known that the class of $k$-trees is an important subclass of trangular graphs. Harry and Plamer [5] first introduced the $2$-tree in 1968, which is showed to be maximal outerplanar graphs in [3, 6]. Beineke and Pippert [1] gave the definition of $k$-trees in 1969. Relating to $k$-trees, there are many interesting applications to the study of a computational complexity and the intersection between graph theory and chemistry [2, 9]. We will just give some notations and definitions below.

\vskip 2mm {\bf Notation 1. } \emph{Let $[a, b]$ be the set of all the integers between $a$ and $b$ with $a \leq b$ including $a, b$, where $a, b$ are integers. Also, let $(a, b] = [a, b] -\{a\}$ and $[a, b) = [a, b] -\{b\}$. In particular, $[a, b] =\phi$ for $a >b$.
}
\vskip 2mm {\bf Notation 2. } \emph{
For ang integer p, if $p \ge 0$, we denote $x_{max\{0, p\}} = x_p$; If $p<0$, we say  $x_{max\{0, p\}} $ does not exist.
}

 \vskip 2mm {\bf Definition 1. } \emph{The $k$-tree, denoted by $T_n^k$, for  positive integers $n, k$ with $n \geq k$, is defined recursively as follows: The smallest $k$-tree is the $k$-clique $K_k$. If $G$ is a $k$-tree with $n \geq k$ vertices and a new vertex $v$ of degree $k$ is added and joined to the vertices of a $k$-clique in $G$, then the obtained graph is a $k$-tree with $n + 1$ vertices.
}

 \vskip 2mm {\bf Definition 2. } \emph{The $k$-path, denoted by $P_n^k$, for positive integers $n, k$ with $n \geq k$, is defined as follows:  Starting with a $k$-clique $G[\{v_1, v_2 ... v_k\}]$. For $i \in [k + 1, n]$, the vetex $v_i$ is adjacent to vertices $\{v_{i - 1}, v_{i - 2} ... v_{i - k} \}$ only.
}

\vskip 2mm {\bf Definition 3. } \emph{The $k$-star, denoted by $S_{k, n-k}$, for positive integers $n, k$ with $n \geq k$, is defined as follows:  Starting with a $k$-clique $G[\{v_1, v_2 ... v_k\}]$ and an independent set $S$ with $|S| = n - k$. For $i \in [k + 1, n]$, the vetex $v_i$ is adjacent to vertices $\{v_{1}, v_{ 2} ... v_{k} \}$ only.
}

\vskip 2mm {\bf Definition 4. } \emph{A vertex $v \in V(T_n^k)$ is called a $k$-simplicial vertex if $v$ is a vertex of degree $k$ whose neighbors form a $k$-clique of $T_n^k$. Let $S_1(T_n^k)$ be the set of all $k$-simplicial vertices of $T_n^k$, for $n \geq k+2$, and set $S_1(K_k) = \phi, S_1(K_{k+1}) = \{v\} $, where $v$ is any vetex of $K_{k+1}$. Let $G =G_0, G_i = G_{i - 1} - v_i$, where $v_i$ is a $k$-simplicial vetex of $G_{i - 1}$, then $\{v_1, v_2 ... v_n \}$ is called a simplicial elimination ordering of the $n$-vertex graph $G$.
}

\vskip 2mm {\bf Definition 5. } \emph{If   $w(G - S) \leq 2$ for   any $k$-clique $G[S]$ of $T_n^k$, we say $T_n^k$ is a hyper pendent edge; If there exists a $k$-clique $G[S]$ with $w(G - S) \geq 3$,   let $C$ be a component of $T_n^k - S$ and contain a unique vertex belonging to $S_1(G)$, then we say that $G[V(S) \cup V(C)]$ is a hyper pendent edge of $T_n^k$, denoted by $\mathcal{P}$. In particular, a $k$-path is a hyper pendent edge.
}

Moreover, let $G[\{v_1, v_2 ... v_k\}]$ denote the initial $k$-clique, then just by the definition of $k$-trees, one can get
\vskip 2mm
{\bf Fact 1. }  \emph  {For the $k$-star, the degree of vertex $v_i$ can be characterized as follows: $d(v_i) = n - k$, for $i \in [1, k]$; $d(v_i) = k$, for $i \in [k+1, n]$. }
\vskip 2mm
{\bf  Fact 2. }  \emph  {For the $k$-path, the degree of vertex $v_i$ can be characterized as follows:\\
 (1) If $4 \leq n \leq 2k$, $d(v_i) = k + i -1 $, for $i \in [1, n - k -1]$; $d(v_i) = n - 1$, for $i \in [n - k, k + 1]$; $d(v_i) = k + n - i$, for $i \in [k + 2, n]$.\\
(2) If $n \geq 2k + 1$, $d(v_i) = k + i - 1$, for $i \in [1, k]$; $d(v_i) = 2k$, for $i \in [k+1, n-k]$; $d(v_i) = k + n -i$, for $i \in [n - k + 1, n]$.}

Easily verified through induction by using the above obseaverations, one can deduce the first generalized multiplicative Zagreb indices and second multiplicative Zagreb indices of the $k$-path and $k$-star as follows.

\vskip 2mm {\bf  Fact 3. }  \emph {
Let $S_{k,n-k}$ be a $k$-star on $n \geq k + 1$ vertices, then \\
$(1)  \prod_{1,c} (S_{k,n-k}) = (n - k)^{ck} k^{c(n-k)};$\\
$(2)  \prod_2 (S_{k,n-k}) = (n - k)^{k(n - k)} k^{k(n - k)}.$}

\vskip 2mm {\bf Fact 4. }  \emph {
Let $P_n^k$ be a $k$-path on $n \geq k + 1$ vertices, then \\
$
(1 . 1)  \prod_{1,c} (P_n^k) = (n - 1)^c \prod_{i = k}^{n - 2} i^{2c} $, if $n \in [k +1, 2k]$;\\
$
(1 . 2)  \prod_{1,c} (P_n^k) = (2k)^{c(n - 2k)} \prod_{i = k}^{2k - 1} i^{2c}$, if $n \geq 2k + 1;$\\
$
(2 . 1)  \prod_2 (P_n^k) = (n - 1)^{n - 1} \prod_{i = k}^{n - 2} i^{2i} $, if $n \in [k +1, 2k]$;\\
$
(2 . 2)  \prod_2 (P_n^k) = (2k)^{2k(n-2k)} \prod_{i = k}^{2k - 1} i^{2i}$, if $n \geq 2k + 1$.}

By considering the derivatives of the following functions, one can get

\vskip 2mm {\bf  Fact 5. }  \emph {
The function $f(x) =\displaystyle{\frac {x}{x+m}}$ is strictly increasing for  $x \in [0, \infty)$, where $m$ is a positive integer.
}

\vskip 2mm {\bf  Fact 6. }  \emph {
The function $f(x) =\displaystyle{\frac {x^x}{(x+m)^{x+m}}}$  is strictly decreasing for  $x \in [0, \infty)$, where $m$ is a positive integer.
}

\section{Main proofs }
Firstly, we give some lemmas that are critical in the proof of our main results.
 \vskip 2mm {\bf Lemma 1 } \emph{    For any $k$-tree $G \not\cong S_{k,n-k}$,  let $u \in S_2$, $N(u) \cap S_1 =\{v_1, v_2... v_s\}$, where $s \geq 1$ is an integer, then\\
(1) For any $i$ with $1 \leq i \leq s$, there exists a vertex $v \in N(u) - \{v_1, v_2 ... v_s\}$ of degree at least $k$ in $G[V(G) - \{v_1, v_2 ... v_s\}]$ such that $vv_i \notin E(G)$.      \\
(2)There exists a $k$-tree $G^*$ such that $\prod_{1,c}(G^*) < \prod_{1,c}(G)$ and  $\prod_2(G^*) > \prod_2(G)$.
  }

{\bf Proof.  } \emph
For (1),  let $G' = G[V(G) - \{v_1, v_2 ... v_s\}]$ and $S = N(u) - \{v_1, v_2 ... v_s\}$,   we obtain that $d_{G'}(u) = |S| = k$ and
$G[S]$ is a $k$-clique by $u \in S_2$.
Since $G \not\cong S_n^k$, $d_{G'}(v) \ge k$ for all $v \in S$.  And by the facts that $N(v_i) \subseteq (N(u) - \{v_1, v_2 ... v_s\}) \cup \{u\}$ with $|N(v_i)| = k$ and $| (N(u) - \{v_1, v_2 ... v_s\}) \cup \{u\}| = k+1$, we have   for any $i \in [1,s]$, there exists a vertex $v \in S$ such that $vv_i \notin E(G)$.

For (2), choose $v_1$ and by (1) there exists a vertex $v \in N(u) - \{v_1, v_2 ... v_s\}$ with $d_{G'}(v) \geq k $ such that $vv_1 \notin E(G)$. If $d_{G'}(v) = k $, and by $uv \in E(G')$, we obtain $G'$ is a $k+1$-clique. Let  $x \in S$ be the vertex such that $d(x) = min_{v\in S} \{d(v)\}$, and let $v_t$ be the vertex such that $v_tx \in E(G)$, $v_ty \notin E(G)$ for some $t \in [1,s]$ and $y \in S$, that is, $ d(x) -1 < d(y)$. Construct a new graph $G^*$ such that $V(G^*) = V(G)$, and $E(G^*) = E(G) - \{v_tx\} + \{v_ty\}$.  Denote $G_0 = G[V(G) - \{x,y\}]$,  since $  d(x) -1 < d(y)$, and by the definition of $\prod_{1,c}(G)$,  $\prod_{2}(G)$ and  Fact 5, Fact 6,  we have

 $$\begin{array}{rcl}
 \displaystyle{ \frac {\prod_{1,c}(G)}{\prod_{1,c}(G^*)}} &=&
 \displaystyle{\frac { [\prod_{w \in V(G_0)} d(w)^c] d(y)^c d(x)^c}{[\prod_{w \in V(G_0)} d(w)^c ] [d(y)+1]^c[d(x)-1]^c}}\\
&=&  \displaystyle{\frac {\displaystyle{d(y)^cd(x)^c}}{\displaystyle{[d(y) + 1]^c[d(x) -1]^c}} } \\
&=&  \displaystyle{\frac { \displaystyle{\frac {d(y)^c}{[d(y) + 1]^c}}} {\displaystyle{\frac{[d(x) - 1]^c} {d(x)^c}} }}\\
&>& 1.
\end{array}$$

Also,

$$\begin{array}{rcl}  \displaystyle{
\frac {\prod_2(G)}{\prod_2(G^*)}}
&=&   \displaystyle{ \frac { [\prod_{w \in V(G_0)}d(w)^{d(w)}] d(y)^{d(y)}d(x)^{d(x)}}{[\prod_{w \in V(G_0)}d(w)^{d(w)} ] [d(y)+1]^{d(y)+1}[d(x)-1]^{d(x)-1}}}  \\
&=&  \displaystyle{ \frac {d(y)^{d(y)}d(x)^{d(x)}}{[d(y)+1]^{d(y)+1}[d(x)-1]^{d(x)-1}}}\\
&=&   \displaystyle{\frac{\displaystyle{[\frac{d(y)^{d(y)}}{[d(y)+1]^{d(y)+1}}]}} { \displaystyle{[\frac {[d(x)-1]^{d(x)-1}}{d(x)^{d(x)}}]}} }\\
&<& 1.\\
\end{array}$$
Thus, we find that the $k$-tree $G^*$ satisfies $\prod_{1,c}(G^*) < \prod_{1,c}(G)$ and  $\prod_2(G^*) > \prod_2(G)$, we are done.

If $d_{G'}(v) \geq k+1 $,  reorder the subindices of $\{v_1, v_2 ... v_s\}$ such that  $vv_i \notin E(G)$ with $i \in [1, s_1]$, where $s_1 \leq s$, and
 by the fact that $G[N(u) - \{v_1, v_2 ... v_s\}]$ is a $k$-clique,  we have $d(u) = k+s$ and $d(v) \geq k+1+s-s_1 $, that is, $ d(v) \geq d(u)  -s_1 +1$.
Construct a new graph $G^*$ such that $V(G^*) = V(G)$, and $E(G^*) = E(G) - \{uv_i\} + \{vv_i\}$, for all $i \in [1,s_1]$.
Since $G[N(u)-\{v_1,v_2 ... v_s\}+ \{u\}]$ is a $k+1$-clique, and for any $i$, $N(v_i) \subseteq N_{G-\{v_1,v_2 ... v_s\}}(u) \cup \{u\}$, then  $G^*$ is a $k$-tree. Denote $G_0 = G[V(G) - \{u, v\}]$,  since $ d(v) \geq d(u)  -s_1 +1$, and by the definition of $\prod_{1,c}(G)$,  $\prod_{2}(G)$ and  Fact 5, Fact 6,  we have

 $$\begin{array}{rcl}
 \displaystyle{ \frac {\prod_{1,c}(G)}{\prod_{1,c}(G^*)}} &=&
 \displaystyle{\frac { [\prod_{w \in V(G_0)} d(w)^c] d(v)^c d(u)^c}{[\prod_{w \in V(G_0)} d(w)^c ] [d(v)+s_1]^c[d(u)-s_1]^c}}\\
&=&  \displaystyle{\frac {\displaystyle{d(v)^cd(u)^c}}{\displaystyle{[d(v) + s_1]^c[d(u) - s_1]^c}} } \\
&=&  \displaystyle{\frac { \displaystyle{[\frac {d(v)^c}{[d(v) + s_1]^c}]}} {\displaystyle{[\frac{[d(u) - s_1]^c} {d(u)^c}]} }}\\
&>& 1.
\end{array}$$

Also,

$$\begin{array}{rcl}  \displaystyle{
\frac {\prod_2(G)}{\prod_2(G^*)}}
&=&   \displaystyle{ \frac { [\prod_{w \in V(G_0)}d(w)^{d(w)}] d(v)^{d(v)}d(u)^{d(u)}}{[\prod_{w \in V(G_0)}d(w)^{d(w)}]  [d(v)+s_1]^{d(v)+s_1}[d(u)-s_1]^{d(u)-s_1}}}  \\
&=&  \displaystyle{ \frac {d(v)^{d(v)}d(u)^{d(u)}}{[d(v)+s_1]^{d(v)+s_1}[d(u)-s_1]^{d(u)-s_1}}}\\
&=&   \displaystyle{\frac{\displaystyle{[\frac{d(v)^{d(v)}}{[d(v)+s_1]^{d(v)+s_1}}]}} { \displaystyle{[\frac {[d(u)-s_1]^{d(u)-s_1}}{d(u)^{d(u)}}]}} }\\
&<& 1.\\
\end{array}$$

Hence, we find that the $k$-tree $G^*$ satisfies $\prod_{1,c}(G^*) < \prod_{1,c}(G)$ and  $\prod_2(G^*) > \prod_2(G)$, we are done.
$\hfill\Box$

\vskip 2mm {\bf Lemma 2 } \emph{Let $G$ be a $k$-tree, if either  $\prod_{1,c}(G)$ attains the maximal or  $\prod_2(G)$ attains the minimal,  then every hyper pendent edge is a k-path.
}

{\bf Proof.  } \emph
Let $\mathcal{P}$$ = G[V(S) \cup V(C)]$ be a hyper pendent edge, where $G[S]=G[\{x_1, x_2 ... x_k\}]$ is a cut $k$-clique and  $V(C) = \{u_1, u_2 ... u_p\}$ with $p$ is a positive ingeter such that  $u_1$ is the only vertex of  $\mathcal{P}$ in $S_1(G)$ and  for $i \in [1, p-1]$, $u_{i}$ is the vertex added following  by $u_{i+1}$ through the process of Definition 1.

\vskip 2mm {\bf  Fact 7. }  \emph {
For any hyper pendent edge  $\mathcal{P}$$ = G[V(S) \cup V(C)]$ as represented above,   $\{u_1, u_2 ... u_{p}\}$ is a simplicial elimination ordering of $\mathcal{P}$.}

{Proof.} By contradiction, assume that  $\{u_1, u_2 ... u_{p}\}$ is not a simplicial elimination ordering of $\mathcal{P}$. Let $u_t$ be the first vertex from $u_1$ to $u_p$ such that $\{u_t, u_{t+1}\} \in S_t$ for $t \in [2, p-1]$, then  $u_{t}u_{t+1} \notin E(G)$ and $\{u_{t}, u_{t+1}\}$ can not be in some $k$-cliques. And by Definition 1, there must be at least two vertices that belongs to $S_1$ in $V(C)$, a contradiction. $\Box$

By Fact 7, we know $\{u_1, u_2 ... u_{p}\}$ is a simplicial elimination ordering of $\mathcal{P}$.  For $p \leq 2$,  $\mathcal{P}$ is a $k$-path by Definition 2; For $p \geq 3$, if  $\mathcal{P}$ is a $k$-path, then we are done. Otherwise, let $u_s$ be the first vertex from $u_p$ to $u_1$ such that  $G[V(S) \cup \{u_{p}, u_{p-1} ... u_{s+1}, u_{s}\}]$ is not a $k$-path. Since $G[V(S) \cup \{u_{p}, u_{p-1} ... u_{s+1}\}]$ is a $k$-path,  for each $i\in[s+1, p]$, let $N_{G-\{u_1, u_2 ... u_{i-1}\}}(u_i) = \{u_{i+1}, u_{i+2} ... ...  u_{min\{p, i+k\}}, x_1, x_2 ... x_{ max\{0, k-p+i\}}\}$, and by Definition 2 and the symmetry of $G[S]$, we have $|N(u_s) \cap \{u_{s+1}, u_{s+2} ... u_{min\{p, s+k\}}\}| = min\{p-s-1,k-1\}$, where $1 \leq s \leq p - 1$.

 For $p \leq k+s$,
suppose that $u_t$ is the vertex such that $u_t \notin N(u_s)$ with $s+2 \leq t \leq p$,  let $N_{G-\{u_1, u_2 ... u_{s-1}\}}(u_s) = \{u_{s+1}, u_{s+2} ... u_{t-1}, u_{t+1} ...  u_{p}, x_1, x_2 ... x_{ k - p + s+1}\} $, and let $|N(x_{k-p+s+1}) \cap \{u_1, u_2 ... u_{s-1}\}| = m $ for $m \in [0, s-1]$. By Definition 2, we have $u_tu_i \notin E(G)$ for $i \in [1, s]$, and then $d(u_{t}) = k + t-s-1$ and $d(x_{k-p+s+1}) > k + p -s +m -1 $. Now construct a new graph $G^*$ such that $V(G^*) = V(G),  E(G^*) = E(G) - \{u_sx_{k-p+s+1}, u_ix_{k-p+s+1}\} + \{u_su_{t}, u_iu_{t}\}$ with $i \in [0, m]$, then $G^*$ is a $k$-tree. Since $t \leq p$, we have $d(x_{k-p+s+1}) > d(u_{t}) +m +1$,  and by the definition of $\prod_{1,c}(G)$, $\prod_2(G)$ and Fact 5, Fact 6,  we get

$$\begin{array}{rcl}
\frac {\prod_{1,c}(G)}{\prod_{1,c}(G^*)} &=& \displaystyle{ \frac {d(u_t)^cd(x_{k-p+s+1})^c}{[d(u_t) +m+1]^c[d(x_{k-p+s+1})-m-1] ^c} } \\
  &=&  \displaystyle{\frac {[\displaystyle{\frac{d(u_t)} {d(u_t)+m+1}}]^c}{[\displaystyle{\frac{d(x_{k-p+s+1}) - m-1}{d(x_{k-p+s+1})}}]^c}} \\
&< & 1,
\end{array}$$
$$\begin{array}{rcl}
\displaystyle{ \frac {\prod_{2}(G)}{\prod_{2}(G^*)}} &=& \displaystyle{ \frac {d(u_t)^{d(u_t)}d(x_{k-p+s+1})^{d(x_{k-p+s+1})}}{[d(u_t) + m+1]^{d(u_t) +m+1}[d(x_{k-p+s+1}) - m-1]^{d(x_{k-p+s+1}) - m-1}} } \\
  &=&  \displaystyle{ \frac {\displaystyle{\frac {d(u_t)^{d(u_t)}}{[d(u_t)+m+1]^{d(u_t)+m+1}}}} {\displaystyle{\frac{[d((x_{k-p+s+1}) - m-1]^{d((x_{k-p+s+1}) - m-1}}{d(x_{k-p+s+1})^{d(x_{k-p+s+1})}}}}} \\
&> & 1.
\end{array}$$

Thus,  $\prod_{1,c}(G^*) > \prod_{1,c}(G)$ and $\prod_{2}(G^*) < \prod_{2}(G)$, a contradiction.

 For $p \geq k+s +1$, let $|N(u_{k+s+1}) \cap \{u_1, u_2 ... u_{s-1}\}| = m$ for $m \in [0, s-1]$.  Since $G[V(S) \cup \{u_{p}, u_{p-1} ... u_{s+1}\}]$ is a $k$-path, we have $G[\{u_{s+1}, u_{s+2} ... u_{s+k+1}\}]$ is a $k+1$-clique. Suppose that $u_t$ is the vertex such that $u_t \notin N(u_s)$ with $s+2 \leq t \leq s+k$, let
 $N_{G- \{u_1, u_2 ... u_{s-1}\}}(u_s) = \{u_{s+1}, u_{s+2} ... u_{t-1}, u_{t+1} ... u_{s+k+1}\}$. Now we construct a new graph  $G^*$ such that $V(G^*) = V(G)$, $E(G^*) = E(G) - \{u_su_{k+s+1}, u_iu_{k+s+1}\} + \{u_su_{t}, u_iu_t\}$ for $i \in [0, m]$,  then $G^*$ is a $k$-tree and $d(u_{k+s+1}) = 2k+ m$, $d(u_{t})= k+t-s-1$.  Since $t \leq s+k$, we have $d(u_{k+s+1}) > d(u_{t})+m+1$, and by the definition of $\prod_{1,c}(G)$, $\prod_2(G)$ and Fact 5, Fact 6,  we get

$$\begin{array}{rcl}
\frac {\prod_{1,c}(G)}{\prod_{1,c}(G^*)} &=& \displaystyle{ \frac {d(u_t)^cd(u_{k+s+1})^c}{[d(u_t) +m+1]^c[d(u_{k+s+1})-m-1] ^c} } \\
  &=&  \displaystyle{\frac{[\displaystyle{\frac {d(u_t)}{d(u_t)+m+1}]^c}} {[\displaystyle{\frac{d(u_{k+s+1})-m - 1}{d(u_{k+s+1})}]^c}}} \\
&< & 1,
\end{array}$$
$$\begin{array}{rcl}
\displaystyle{ \frac {\prod_{2}(G)}{\prod_{2}(G^*)}} &=& \displaystyle{ \frac {d(u_t)^{d(u_t)}d(u_{k+s+1})^{d(u_{k+s+1})}}{[d(u_t) +m+ 1]^{d(u_t)+m +1}[d(u_{k+s+1})-m - 1]^{d(u_{k+s+1})-m - 1}} } \\
  &=&  \displaystyle{ \frac {\displaystyle{\frac{d(u_t)^{d(u_t)}} {[d(u_t)+m+1]^{d(u_t)+m+1}}}}{\displaystyle{\frac{[d((u_{k+s+1})-m - 1]^{d(u_{k+s+1}) -m- 1}}{d(u_{k+s+1})^{d(u_{k+s+1})}}}}} \\
&> & 1.
\end{array}$$

Thus,  $\prod_{1,c}(G^*) > \prod_{1,c}(G)$ and $\prod_{2}(G^*) < \prod_{2}(G)$, a contradiction.
Hence, for any $s \in [1, p]$ $N_{G-\{u_1, u_2 ... u_{s-1}\}}(u_s) = \{u_{s+1}, u_{s+2} ... u_{min \{p, k+s\}}, x_1, x_2 ... x_{max \{0, k - p + s\}}\} $, that is,  $\mathcal{P}$ is a $k$-tree. $\hfill\Box$

\vskip 2mm {\bf Lemma 3 } \emph{    Let $G$ be a  $k$-tree, if either $\prod_{1,c}(G)$ attains the maximal or  $\prod_2(G)$ attains the minimal, then $|S_1(G)| = 2$.
 }

{\bf Proof.  } \emph
We know that $|S_1(G)| \geq 2$ for $n \geq k + 2$, and by Lemma 2, every hyper pendent edge is a $k$-path  for $\prod_{1,c}(G)$ to attain the maximal or  $\prod_2(G)$ to attain the minimal. If $|S_1(G)| = 2$, we are done; Suppose that $|S_1(G)| \geq 3$, it suffices to prove that there exists a graph $G'$ such that $|S_1(G')| = |S_1(G)| - 1$ with $ \prod_{1,c}(G') > \prod_{1,c}(G)$ and  $\prod_2(G') < \prod_2(G)$.

\vskip 2mm {\bf  Fact 8. }  \emph {
 For any $k$-tree $G$ satisfying the conditions of Lemma 3, if $|S_1(G)| \geq 3$, then there exists a $k$-clique $G[S]$ such that $w(G-S) \geq 3$.}

{Proof.}  We will proceed by induction on $n = |G|$. For $n= k+3$, it's trivial; For $n \geq k+4$, assume that the fact is true for the $k$-tree $G$ with  $n < k+p$, and consider $n = k+p$.
 If $|S_1(G)| \ge 4$, choose any vertex $v \in S_1(G)$, or $|S_1(G)| = 3$ and $|S_2(G)| \ge 2$, choose the vertex $v \in S_1(G)$ such that $ N(w) \cap S_1(G) = \{v\}$ for some $w \in S_2(G)$, then construct a new graph $G'$ such that $G' = G - v$. Since $S_2(G)$ is an dependent set and $G[N(v)]$ is a $k$-clique for any $v \in S_1(G)$, we obtain $|S_1(G')| \ge 3$. By the induction hypothesis, there exists a $k$-clique $G[S]$ in $G'$ such that $w(G' - S) \ge 3$.
Thus,  by adding back $v$, $G[S]$ is still a $k$-clique in $G$ and $w(G-S) \geq 3$, we are done. Next, we only consider  $|S_1(G)| = 3$ and $|S_2(G)| = 1$.

 Let $S_1(G) = \{v_1, v_2, v_3\}$ and $G_0 = G-\{v_1, v_2, v_3\}$,  by Definition 4, we have $G_0$  is a $k+1$-clique, denoted $G[\{x_1, x_2 ... x_{k+1}\}]$. If there exists $N(v_i) = N(v_j)$,  for some $i, j \in [1, 3]$ with $i\neq j$, and take $S = N(v_i)$, then $w(G- S) \ge 3$, we are done;  If $N(v_i) \neq N(v_j)$,  for any $i, j \in [1, 3]$ with $i \neq j$, then reorder the index of $x_i$ such that $N(v_1) = \{x_1, x_2 ... x_k\}$, $N(v_2) = \{x_2, x_3 ... x_{k+1}\}$ and $N(v_3) = \{x_1, x_3 ... x_{k+1}\}$. Construct a new graph $G^*$ such that $V(G^*) = V(G)$, $E(G^*) = E(G) - \{v_1x_1\} +\{v_1v_2\}$, then $G^*$ is still a $k$-tree and $d_G(x_1) = k+2$,  $d_{G^*}(x_1) = k+1$, $d_G(v_1) =d_G(v_2) =k$ and  $d_{G^*}(v_2) = k+1$.  By the definition of $\prod_{1,c}(G)$, $\prod_2(G)$ and Fact 6, we have
$$\begin{array}{rcl}
\displaystyle{\frac {\prod_{1,c}(G)}{\prod_{1,c}(G^*)}} &=& \displaystyle{ \frac {d(v_2)^cd(x_{1})^c}{[d(v_2) +1 ]^c[d(x_{1})-1] ^c} } \\
  &=&  \displaystyle{[\frac {k(k+2)}{(k+1)^2}]^c} \\
&< & 1,
\end{array}$$
$$\begin{array}{rcl}
\displaystyle{ \frac {\prod_{2}(G)}{\prod_{2}(G^*)}} &=& \displaystyle{ \frac {d(v_2)^{d(v_2)}d(x_{1})^{d(x_{1})}}{[d(v_2) +1]^{d(v_2)+1 }[d(x_{1})-1]^{d(x_{1})-1 }} } \\
  &=&  \displaystyle{ \frac {(k+2)^{k+2}k^k}{(k+1)^{2(k+1)}}} \\
&=& \displaystyle{ \frac {\displaystyle{[\frac {k^k}{(k+1)^{k+1}}]} }{\displaystyle{[\frac {(k+1)^{k+1}}{(k+2)^{k+2}}]} }  } \\
&> & 1.
\end{array}$$

Thus, we find a graph $G^*$ with $\prod_{1,c}(G^*) > \prod_{1,c}(G)$ and $\prod_{2}(G^*) < \prod_{2}(G)$, a contradiction with that $\prod_{1,c}(G)$ attains the maximal or  $\prod_2(G)$ attains the minimal, we are done.
   $\Box$

Choose a $k$-clique $G[S]$ with $w(G - S) \geq 3$ such that there are two components of $G - S$: $C_1, C_2$ with $|C_1|=p, |C_2| = q$ and $p + q$ being minimal, for $p \geq q$. Let $u_1\in V(C_1)$,  $v_1\in V(C_2)$ with $\{u_1, v_1\} \subseteq S_1(G)$. Let $N_{G- \{v_1,v_2 ... v_{i-1}\}}(v_i) = \{v_{i+1}, v_{i+2} ... v_{min\{k+1, q\}}, x_1, x_2 ... x_{max\{0, {k-q+i}\}}$,
$N_{G-\{u_1,u_2 ... u_{j-1}\}}(u_j) = \{u_{j+1}, u_{j+2} ... u_{min\{k+1, p\}}, y_1,y_2 ... y_{max\{0, {k-p+i}\}}\}$ for $i \geq 1, j \geq 1$, where $\{v_1, v_2 ... v_q\}$ and $\{u_1, u_2 ... u_p\}$ are simplicial elimination orderings of $G[S\cup V(C_1)]$ and $G[S\cup V(C_2)]$, respectively.  We will prove Lemma 3 by induction on $q$.

$(1)$ If $q =1$, then  $d(v_1) = k$. Choose $x_t \in N(v_1)$, let $|N(x_t) \cap \{u_1, u_2 ... u_p\}| = m$ for $m \in [1, k]$,  we get $d(x_t) > k+1+m$ by $w(G -S) \geq 3$, and then $d(x_t) > d(v_1) + m +1$. Now construct a new graph  $G^*$ such that $V(G^*) = V(G), E(G^*) = E(G) - \{u_ix_t\} + \{u_iv_{1}\}$ for $i \in [1, m]$,  then $G^*$ is a $k$-tree and $|C_1| + |C_2| = p$ with $G[\{x_1, x_2 ... x_{t-1}, x_{t+1} ... x_k, v_1\}]$ is a $k$-clique in $G^*$.  Since   $d(x_{t}) > d(v_{1})+m+1$,  by the definition of $\prod_{1,c}(G)$, $\prod_2(G)$ and Fact 5, Fact 6, we have

$$\begin{array}{rcl}
\frac {\prod_{1,c}(G)}{\prod_{1,c}(G^*)} &=& \displaystyle{ \frac {d(v_1)^cd(x_{t})^c}{[d(v_1) +m ]^c[d(x_{t})-m] ^c} } \\
  &=&  \displaystyle{\frac {[\displaystyle{\frac {d(v_1)}{d(v_1)+m}]^c}}{[\displaystyle{\frac{d(x_{t})-m }{d(x_{t})}]^c}}} \\
&< & 1,
\end{array}$$
$$\begin{array}{rcl}
\displaystyle{ \frac {\prod_{2}(G)}{\prod_{2}(G^*)}} &=& \displaystyle{ \frac {d(v_1)^{d(v_1)}d(x_{t})^{d(x_{t})}}{[d(v_1) +m]^{d(v_1)+m }[d(x_{t})-m]^{d(x_{t})-m }} } \\
  &=&  \displaystyle{ \frac {\displaystyle{\frac{d(v_1)^{d(v_1)}} {[d(v_1)+m ]^{d(v_1)+m }}}}{\displaystyle{\frac{[d((x_{t})-m ]^{d((x_{t}) -m}}{d(x_{t})^{d(x_{t})}}}}} \\
&> & 1.
\end{array}$$

Then,  $\prod_{1,c}(G^*) > \prod_{1,c}(G)$ and $\prod_{2}(G^*) < \prod_{2}(G)$. Thus,  let $G' = G^*$,  $|S_1(G')| = |S_1(G)| - 1$, $ \prod_{1,c}(G') > \prod_{1,c}(G)$ and  $\prod_2(G') < \prod_2(G)$, and we are done.

$(2)$ Assume that $q =s $, there exists a $k$-tree $G'$ such that $|S_1(G')| = |S_1(G)| -1$, $ \prod_{1,c}(G') > \prod_{1,c}(G)$,  $\prod_2(G') < \prod_2(G)$ and we consider $q = s+1 $.

If  $q \leq k$,  we have  $d(v_q) = k+q-1$ by the fact that $G[S\cup V(C_2)]$ is a $k$-path.
 Choose $x_t \in N(v_1)$,  we know $x_t \in N(v_i)$ for all $i \in [1, p]$ by  $G[S\cup V(C_2)]$ is a $k$-path. Let $|N(x_t) \cap \{u_1, u_2 ... u_p\}| = m$ for $m \in [1, k]$,  we have  $d(x_t) > k+q+m$ by $w(G - S) \geq 3$, and then $d(x_t) > d(v_q) + m +1$. Now construct a new graph  $G^*$ such that $V(G^*) = V(G), E(G^*) = E(G) - \{u_ix_t\} + \{u_iv_{q}\}$ for $i \in [1, m]$,  then $G^*$ is a $k$-tree and $|C_1| + |C_2| = p+q-1$ with $G[\{x_1, x_2 ... x_{t-1}, x_{t+1} ... x_k, v_q\}]$ is a $k$-clique in $G^*$.  Since   $d(x_{t}) > d(v_{q})+m+1$, by the definition of $\prod_{1,c}(G)$, $\prod_2(G)$ and Fact 5, Fact 6, we have

$$\begin{array}{rcl}
\frac {\prod_{1,c}(G)}{\prod_{1,c}(G^*)} &=& \displaystyle{ \frac {d(v_q)^cd(x_{t})^c}{[d(v_q) +m ]^c[d(x_{t})-m] ^c} } \\
  &=&  \displaystyle{\frac{[\displaystyle{\frac {d(v_q)}{d(v_q)+m}}]^c} {[\displaystyle{\frac{d(x_{t})-m }{d(x_{t})}}]^c}} \\
&< & 1,
\end{array}$$
$$\begin{array}{rcl}
\displaystyle{ \frac {\prod_{2}(G)}{\prod_{2}(G^*)}} &=& \displaystyle{ \frac {d(v_q)^{d(v_q)}d(x_{t})^{d(x_{t})}}{[d(v_q) +m]^{d(v_q)+m }[d(x_{t})-m]^{d(x_{t})-m }} } \\
  &=&  \displaystyle{ \frac {\displaystyle{\frac {d(v_q)^{d(v_q)}}{[d(v_q)+m ]^{d(v_q)+m }}}}{\displaystyle{\frac{[d((x_{t})-m ]^{d((x_{t}) -m}}{d(x_{t})^{d(x_{t})}}}}} \\
&> & 1.
\end{array}$$

Then, $ \prod_{1,c}(G) < \prod_{1,c}(G^*)$, $\prod_2(G) > \prod_2(G^*)$ and $q =s$ in $G^*$, then by the induction hypothesis, there exists a $k$-tree $G'$ such that $|S_1(G')| = |S_1(G)| -1$, $ \prod_{1,c}(G') > \prod_{1,c}(G)$ and  $\prod_2(G') < \prod_2(G)$, we are done.

 If $q \geq k+1$,   we have  $N(u_1) = \{u_2, u_3 ... u_{k+1}\}$, $N(v_1) = \{v_2, v_3 ... v_{k+1}\}$ by the facts that $p\geq q$ and $G[S\cup V(C_1)]$, $G[S\cup V(C_2)]$ are  $k$-paths. We construct a new graph $G^*$ such that  $V(G^*) = V(G), E(G^*) = E(G) - \{v_1v_i\} + \{u_jv_{1}\}$ for $i \in [2, k+1]$, $j \in [1, k]$. And by Fact 2 and the definition of $\prod_{1,c}(G)$, $\prod_2(G)$,
we obtain
$$\begin{array}{rcl}$$
 \displaystyle{\frac {\prod_{1,c}(G)}{\prod_{1,c}(G^*)}} &=&
\displaystyle{ \frac { \prod_{i=2}^{k+1}d(v_i)^c \prod_{j=1}^{k}d(u_j)^c}{ \prod_{i=2}^{k+1}[d(v_i)-1]^c
\prod_{j=1}^{k, }[d(u_j)+1]^c}}     \\
&=& 1,\\
 \displaystyle{\frac {\prod_2(G)}{\prod_2(G^*)}} &=&
\displaystyle{  \frac { \prod_{i=2}^{k+1}d(v_i)^{d(v_i)} \prod_{j=1}^{k}d(u_j)^{d(u_j)}}{ \prod_{i=2}^{k+1}[d(v_i)-1]^{d(v_i)-1}
\prod_{j=1}^{k}[d(u_j)+1]^{d(u_j)+1}}
}\\
&=& 1.
\end{array}$$

Then, $ \prod_{1,c}(G) = \prod_{1,c}(G^*)$, $\prod_2(G) = \prod_2(G^*)$ and $q =s$ in $G^*$, then by the induction hypothesis, there exists a $k$-tree $G'$ such that $|S_1(G')| = |S_1(G)| -1$, $ \prod_{1,c}(G') > \prod_{1,c}(G)$ and  $\prod_2(G') < \prod_2(G)$, we are done.
 $\hfill\Box$

Now, we turn to prove the main results of the paper.

{\bf Proof of Theorem 2.  } \emph
For any $k$-tree $T_n^k$, if $|S_1(T_n^k)| = n - k$, then $T_n^k \cong S_{k,n-k}$, we are done.
And if $|S_1(T_n^k)| \leq n - k -1$, we can recursively use Lemma 1 to make $\prod_{1,c}(T_n^k)$ decreasing until $|S_1(T_n^k)| = n - k$. Thus, we have  $T_n^k \cong S_{k,n-k}$ for  $\prod_{1,c}(T_n^k)$ to arrive the minimal value.

By Lemma 2, if $\prod_{1,c}(T_n^k)$ get the maximal, then every hyper pendent edge is a $k$-path, and by Lemma 3, $|S_1(T_n^k)| = 2$, implying that  $T_n^k \cong P_n^k$ for  $\prod_{1,c}(T_n^k)$ to arrive the maximal value.    $\hfill\Box$

 \vskip 2mm
{\bf Proof of Theorem 3.  } \emph
For any $k$-tree $T_n^k$, if $|S_1(T_n^k)| = n - k$, then $T_n^k \cong S_{k,n-k}$ ,we are done.
And if $|S_1(T_n^k)| \leq n - k -1$, we can recursively use Lemma 1 to make $\prod_2(T_n^k)$ increasing until $|S_1(T_n^k)| = n - k$, then we have $T_n^k \cong S_{k,n-k}$ for  $\prod_2(T_n^k)$ to arrive the maximal value.

By Lemma 2, if $\prod_2(T_n^k)$ get the minimal,  every hyper pendent edge is a $k$-path, and by Lemma 3, $|S_1(T_n^k)| = 2$. Then this $k$-tree is a $k$-path, that is, $T_n^k \cong P_n^k$ for  $\prod_2(T_n^k)$ to arrive the minimal value.
 $\hfill\Box$


\section*{The Methods}

Our Monte Carlo evolutionary modeling is performed with the {\tt StarTrack} binary
population synthesis code~\cite{Belczynski2002}. In particular, we incorporate a
calibrated treatment of tidal interactions in close binaries~\cite{Belczynski2008a},
a physical measure of the common envelope (CE) binding energy~\cite{Dominik2012,Xu2010},
and a rapid explosion supernova model that reproduces the observed mass gap between
neutron stars and black holes~\cite{Fryer2012,Belczynski2012}. Our updated mass spectrum
of black holes shows a strong dependence on the metallicity of the progenitor stars
(Extended Data Fig.~\ref{fig:BHMassRel}). In galaxies with metallicities similar to the
Milky Way ($Z={\rm ~Z}_{\odot}=0.02$) black holes formed out of single massive stars (initial mass
$M_{\rm ZAMS}=150{\rm ~M}_{\odot}$) reach a maximum mass of $M_{\rm BH}=15{\rm ~M}_{\odot}$, while for very
low metallicity ($Z=0.0001=0.5\%{\rm ~Z}_{\odot}$) the maximum mass becomes $M_{\rm BH}=94{\rm ~M}_{\odot}$.
The above input physics represents our standard model (M1), which is representative of
our classical formation scheme for double compact objects (BH-BH, BH-NS, and NS-NS).

We have adopted specific values for a number of evolutionary
parameters. Single stars are evolved with calibrated formulae based on detailed
evolutionary calculations~\cite{Hurley2000}. Massive star winds are adopted from
detailed studies of radiation driven mass loss~\cite{Vink2011}. For the Luminous Blue
Variable phase the high mass loss rate is adopted ($1.5 \times 10^{-4} {\rm ~M}_{\odot}$ yr$^{-1}$).
Binary interactions, and in particular the stability of RLOF, is judged based on
binary parameters: mass ratio, evolutionary stage of donor, response to mass loss, and
behavior of the orbital separation in response to mass transfer. The orbital separation
is additionally affected by gravitational radiation, magnetic braking, and angular
momentum loss associated with systemic mass loss. During stable RLOF we assume that
half of the mass is accreted onto the companion, while the other half ($1-f_{\rm a}=0.5$)
is lost with specific angular momentum ($dJ/dt= j_{\rm loss}
[J_{\rm orb}/(M_{\rm don}+M_{\rm acc})] (1-f_{\rm a}) dM_{\rm RLOF}/dt$ with
$j_{\rm loss}=1.0$~\cite{Podsiadlowski1992}).
CE is treated with energy balance with fully effective conversion of orbital energy into
envelope ejection ($\alpha=1.0$), while the envelope binding energy for massive stars
is calibrated by a parameter $\lambda$ that depends on star radius, mass, and metallicity.
For massive stars $\lambda \approx 0.1$ is adopted~\cite{Xu2010}.
During CE compact objects accrete at $10\%$ Bondi-Hoyle rate as estimated by recent
hydrodynamical simulations~\cite{Ricker2008,Macleod2015}.
Our CE is done ``instantaneously'', so the time at the beginning
and end of CE is exactly the same (see Fig.~\ref{fig:ExamplEvol}); the
time duration of CE has no impact on our results.

We consider two extra variations of the binary evolution input physics. In one model
(M2) we test highly uncertain CE physics~\cite{Ivanova2013} and we allow for
Hertzsprung gap (HG) stars to initiate and survive CE evolution. This is an
optimistic assumption, since these stars may not allow for CE
evolution~\cite{Pavlovskii2015}, nor survive as a binary if CE does
happen~\cite{Belczynski2007}.
For comparison, in our standard model we allow only evolved stars with a deep
convective envelope (core Helium burning stars) to survive CE.

In the opposite extreme, we employ a model (M3) where black holes receive high
natal kicks. In particular, each BH gets a natal kick with its
components drawn from a $1$-D Maxwellian distribution with $\sigma=265$ km
s$^{-1}$, independent of BH mass.
Such high natal kicks are measured for Galactic pulsars~\cite{Hobbs2005}.
This is a pessimistic assumption, as high natal kicks tend to disrupt BH-BH
progenitor binaries. This assumption is not yet excluded based on
electromagnetic observations~\cite{Belczynski2015}.
In contrast, in our standard model BH natal kicks decrease with BH mass. In
particular, for massive BHs that form through direct collapse of an entire
star to a BH with no supernova explosion
($M_{\rm BH} \gtrsim 10{\rm ~M}_{\odot}$ for solar metallicity;
$M_{\rm BH} \gtrsim 15{\rm ~M}_{\odot}$ for $Z=10\%{\rm ~Z}_{\odot}$; and $M_{\rm BH} \gtrsim 15$--$30{\rm ~M}_{\odot}$
for $Z=1\%{\rm ~Z}_{\odot}$) we assume no natal kicks~\cite{Fryer2012}.
We have also calculated a series of models with intermediate BH kicks (see Extended Data
Fig.~\ref{fig:IntrRate}):
\mbox{$\sigma=200$ km\ s}$^{-1}$ (model M4),
$\sigma=130$ km\ s$^{-1}$ (model M5),
$\sigma=70$ km\ s$^{-1}$ (model M6).

For each evolutionary model we compute $2 \times 10^7$ massive binaries for each
point on a grid of $32$ sub-models covering a wide range of metallicities:
$Z=0.0001$, $0.0002$, $0.0003$, $0.0004$, $0.0005$, $0.0006$, $0.0007$, $0.0008$,
$0.0009$, $0.001$, $0.0015$, $0.002$, $0.0025$, $0.003$, $0.0035$, $0.004$, $0.0045$,
$0.005$, $0.0055$, $0.006$, $0.0065$, $0.007$, $0.0075$, $0.008$, $0.0085$, $0.009$,
$0.0095$, $0.01$, $0.015$, $0.02$, $0.025$, $0.03$.
We assume that stellar evolution at even lower metallicities proceeds in the same way
as the evolution at $Z=0.5\%{\rm ~Z}_{\odot}$. However, we note that stars with very low metal
content (e.g., Population III) may evolve differently than metal-rich
stars~\cite{Szecsi2015}.

Each sub-model is computed with initial distributions of orbital periods
($\propto (\log P)^{-0.5}$), eccentricities ($\propto e^{-0.42}$), and mass ratios
($\propto q^{0}$) appropriate for massive stars~\cite{Sana2012}.
We adopt an initial mass function that is close to flat for low mass stars
($\propto M^{-1.3}$ for $0.08 \leq M<0.5{\rm ~M}_{\odot}$ and $\propto M^{-2.2}$ for $0.5
\leq M<1.0{\rm ~M}_{\odot}$) and top heavy for massive stars ($\propto M^{-2.3}$ for $1.0 \leq
M \leq 150{\rm ~M}_{\odot}$), as guided by recent observations~\cite{Bastian2010}. The adopted
IMF generates higher BH-BH merger rate densities as compared with the steeper IMF
($\propto M^{-2.7}$ for $1.0 \leq M \leq 150{\rm ~M}_{\odot}$) adopted in our earlier
studies~\cite{Dominik2013,Belczynski2015} as there are more BH-BH merger
progenitors in our simulations~\cite{deMink2015}.

A moderate binary fraction ($f_{\rm bi}=0.5$) is adopted for stars with mass
$M_{\rm ZAMS}<10{\rm ~M}_{\odot}$, while we assume that all more massive stars are formed in
binaries ($f_{\rm bi}=1.0$) as indicated by recent empirical
estimates~\cite{Sana2012,Duchene2013}).

We adopt an extinction corrected cosmic star formation rate based on numerous
multi-wavelength observations~\cite{Madau2014}:
\begin{equation}
\mbox{SFR}(z)=0.015 {(1+z)^{2.7} \over 1+[(1+z)/2.9]^{5.6}} \,{\rm ~M}_{\odot} \, {\rm Mpc}^{-3}
\, {\rm yr}^{-1}.
\label{sfr}
\end{equation}
This SFR declines rapidly at high redshifts ($z>2$). This may be contrasted with some
SFR models that we have used in the past~\cite{Strolger2004} which generated a greater
number of stars at high redshifts. This revision will thus reduce the BH-BH merger
rate densities at {\em all} redshifts. Even though the formation of BH-BH binaries
takes a very short time ($\sim 5$ Myr), the time to coalescence of two black holes
may take a very long time (Fig.~\ref{fig:ExamplEvol} and Extended Data
Fig.~\ref{fig:Emergence}).

In our new treatment of chemical enrichment of the Universe we follow the mean
metallicity increase with cosmic time (since Big Bang till present). The mean
metallicity as a function of redshift is given by
\begin{equation}
\log(Z_{\rm mean}(z))=0.5+\log \left(  {y \, (1-R) \over \rho_{\rm b}}
\int_z^{20} {97.8\times10^{10} \, sfr(z') \over H_0 \, E(z') \, (1+z')} dz' \right)
\label{Zmean}
\end{equation}
with a return fraction $R=0.27$ (mass fraction of each generation of stars that
is put back into the interstellar medium), a net metal yield $y=0.019$ (mass of new
metal created and ejected into the interstellar medium by each generation of stars
per unit mass locked in stars), a baryon density $\rho_{\rm b}=2.77 \times 10^{11}
\,\Omega_{\rm b}\,h_0^2\,{\rm ~M}_{\odot}\,{\rm Mpc}^{-3}$ with $\Omega_{\rm b}=0.045$
and $h_0=0.7$, a star formation rate given by eq.~\ref{sfr}, and
$E(z)=\sqrt{\Omega_{\rm M}(1+z)^3+\Omega_{\rm k}(1+z)^2+\Omega_\Lambda)}$
with $\Omega_\Lambda=0.7$, $\Omega_{\rm M}=0.3$, $\Omega_{\rm k}=0$,
and $H_0=70.0\,{\rm km}\,{\rm s}^{-1}\,{\rm Mpc}^{-1}$.
The shape of the mean metallicity dependence on redshift follows recent
estimates~\cite{Madau2014}, although the level was increased by $0.5$ dex to better
fit observational data~\cite{Vangioni2015}.
At each redshift we assume a log--normal distribution of metallicity around the
mean, with $\sigma=0.5$ dex~\cite{Dvorkin2015}.
Our prescription (Extended Data Fig.~\ref{fig:MetEvol}) produces more low-metallicity
stars than previously~\cite{Dominik2013}. Since BH-BH formation is enhanced at
low-metallicity~\cite{Belczynski2010a}, our new approach increases the predicted rate
densities of BH-BH mergers.

Here we discuss caveats of evolutionary calculations.
First, we only consider isolated binary evolution, and thus our approach is
applicable to field stars in low density environments. It is possible that
dynamical interactions enhance BH-BH merger formation in dense globular
clusters~\cite{Rodriguez2016}, offering a completely independent channel.

Second, our predictions are based on a ``classical'' theory of stellar and
binary evolution for the modeling of massive stars that we have compiled,
developed, and calibrated over the last $15$ years. We do not consider exotic
channels for the formation of BH-BH mergers, such as the one from rapidly
rotating stars in contact binaries~\cite{Almeida2015}.

Third, our modeling includes only three evolutionary models: a ``standard''
model consisting of our best estimates for reasonable parameters, as well as
``optimistic'' and ``pessimistic'' alternate models. The optimistic model
consists of only one change from the standard model: we allow all stars beyond
the main sequence to survive the common envelope phase. Alternatively, the
pessimistic model also consists of only one change: larger BH natal kicks.
We have not investigated other possible deviations from the standard model (e.g.,
different assumptions of mass and angular momentum loss during stable mass
transfer evolution) nor have we checked inter-parameter degeneracies (e.g.,
models with high BH kicks {\em and}\/ an optimistic common envelope phase).
Albeit with low statistics and limited scope, precursor versions of these
computationally demanding studies have already been performed~\cite{O'Shaughnessy2007};
these calculations indicate that our three models are likely to cover the range
of interesting effects.

Fourth, our observations are severely statistically limited. We are attempting
to draw inferences about our models based on a single detection (GW150914).

In was argued~\cite{AstroPaper} that the formation of GW150914 in isolated
binary evolution requires a metallicity lower than $50\%{\rm ~Z}_{\odot}$. This was based on
single stellar models~\cite{Spera2015}; stars in close binaries are subject to
significant mass loss during RLOF/CE, and they form BHs with lower mass than single
stars. Thus in binaries the metallicity threshold for massive BH formation is lower
than in single stellar evolution. For example, formation of a {\em single}\/ $30{\rm ~M}_{\odot}$
BH requires $Z<25\%{\rm ~Z}_{\odot}$ (Extended Data Fig.~\ref{fig:BHMassRel}), while formation
of two such BHs in a {\em binary} requires $Z<10\%{\rm ~Z}_{\odot}$ (Extended Data
Fig.~\ref{fig:BHBHMass}). The value of this threshold depends on assumptions
for the model of stellar evolution, winds, and BH formation processes. The physical
models we have adopted yield a threshold of $Z<10\%{\rm ~Z}_{\odot}$, the same as the one
obtained with MESA for homogeneous stellar evolution~\cite{Marchant2016}.
Our model was calibrated on known masses of BHs, and in particular we do not
exceed $15{\rm ~M}_{\odot}$ for ${\rm ~Z}_{\odot}$ (the highest mass stellar BH known in our Galaxy).
In contrast, single stellar models used to derive the high metallicity threshold
produce $25{\rm ~M}_{\odot}$ for ${\rm ~Z}_{\odot}$~\cite{Spera2015}. The highest threshold obtained
with binary evolution was reported at the level of $50\%{\rm ~Z}_{\odot}$~\cite{Eldridge2016}.
Such high value of metallicity threshold for the progenitor of GW150914 implies that
stars at approximately solar metallicity ($Z=0.014$) produce BHs as massive as
$40{\rm ~M}_{\odot}$~\cite{Eldridge2016}. This is neither supported nor excluded by available
electromagnetic BH mass measurements~\cite{fn2}.

In the following we present calculation of the gravitational radiation signal.
The output of {\tt StarTrack} is a binary merger at a given time. We then
calculate the gravitational waveform associated with this merger, and determine
whether this binary would have been observable by LIGO in the O1
configuration~\cite{Dominik2015,Belczynski2015}.

We model the full inspiral-merger-ringdown waveform of the binaries using the
IMRPhenomD gravitational waveform template family~\cite{Khan2015,Husa2015}. This
is a simple and fast waveform family which neglects the effects of spin (which
are not relevant for GW150914). We consider a detection to be given by a threshold
SNR $>8$ in a single detector, and we use the fiducial O1 noise curve~\cite{fn1}.
We calculate the face-on, overhead SNR for each binary directly from Eq.~2 of
\cite{Dominik2015}. We then calculate the luminosity
distance at which this binary would be detected with $\mbox{SNR}=8$. Note that
as the distance to the binary changes, the observer frame (redshifted) mass also
changes, and therefore calculating the horizon redshift requires an iterative
process. Once this has been calculated, we then determine the predicted
detection rates using Eq.~9 of~\cite{Dominik2015}. The effects of the antenna
power pattern are incorporated in the $p_{\rm det}$ term in this equation.

Estimate of fiducial aLIGO sensitivity during the 16-day GW150914 analysis
is shown in Figure~\ref{fig:RateComp}. We estimate the sensitivity to coalescing
compact binaries using a reference O1 noise curve. We assume that both detectors
operate with the fiducial O1 noise curve, which is the same sensitivity we adopted
to calculate compact binary detection rates. For comparison, this model agrees
reasonably well with the ``early-high'' sensitivity model~\cite{LVC2013}. Our
expression is a 50th percentile upper limit, assuming no detections. The critical
application of this expression is not related to its overall normalization; we
are instead interested in its shape, which characterizes the strongly mass-dependent
selection biases of LIGO searches.

Using these inputs, our fiducial estimate of the advanced LIGO sensitivity during
the first 16 days of O1 for a specific mass bin,
$\Delta M_{\rm i}$, is given by
\begin{eqnarray}
R_{{\rm D},\Delta M_{\rm i}, UL} = \frac{0.7}{V_{\Delta M_{\rm i}}T} ,
\end{eqnarray}
where $T$ is 16 days, corresponding to the analysis of GW150914~\cite{DiscoveryPaper},
and the volume
$$
V_{\Delta M_{\rm i}} = (\Delta M_{\rm i})^{-1}\int dM \int \frac{dz}{1+z} \frac{dV}{dz}
p_{\rm det}(<w,M)
$$
is the  sensitive volume averaged over the mass bin $\Delta M_{\rm i}$, and
$p_{\rm det}(<w,M)$ is the orientation-averaged detection
probability~\cite{Dominik2015,Belczynski2015}.
The function $p_{\rm det}(<w,M)$ depends on the coalescing binary redshifted mass
through the maximum luminosity distance (``horizon distance'') at which a source could
produce a response of SNR$>$8 in a single detector. To calculate this distance, we
adopt the same model for inspiral, merger, and ringdown~\cite{Husa2015,Khan2015} used
in the text to estimate compact binary detection rates. Extended Data
Figure~\ref{fig:Horizon} shows our estimated horizon redshift as a function of the
total redshifted binary merger mass for equal mass mergers.

\noindent
{\bf Code availability.}
We have opted not to release the population synthesis code {\tt StarTrack}
used to generate binary populations for this study.