# On Correlation Coefficients

Alexei Stepanov Education and Research Cluster “Institute of High Technology”,
Immanuel Kant Baltic Federal University, A.Nevskogo 14, Kaliningrad, 236041
Russia, email<EMAIL_ADDRESS>
###### Abstract

In the present paper, we discuss the Pearson $\rho$, Spearman $\rho_{S}$,
Kendall $\tau$ correlation coefficients and their statistical analogues
$\rho_{n},\rho_{n,S}$ and $\tau_{n}$. We propose a new correlation coefficient
$r$ and its statistical analogue $r_{n}$. The coefficient $r$ is based on
Kendal’s and Spearman’s correlation coefficients. A new extension of the
Pearson correlation coefficient is also discussed. We conduct simulation
experiments and study the behavior of the above correlation coefficients. We
observe that the behavior of $\rho_{n}$ can be very different from the
behavior of the rank correlation coefficients $\rho_{n,S},\tau_{n}$ and
$r_{n}$, which, in turn, behave in a similar way. The question arises: which
correlation coefficient better measures the dependence rate? We try to answer
this question in the final conclusion.

Keywords and Phrases: bivariate distributions; Pearson, Kendall and Spearman
correlation coefficients.

AMS 2000 Subject Classification: 60G70, 62G30

##  1 Introduction

Let $(X,Y),(X_{1},Y_{1}),\ldots,(X_{n},Y_{n})$ be independent and identically
distributed random vectors with absolutely continuous bivariate distribution
function $F(x,y)$, density function $f(x,y)$ and marginal distribution
functions $H(x)$ and $G(y)$. A rate of dependence between $X$ and $Y$
indicates how closely these variables are related. This ranges from complete
dependence to independence. There are three basic measures of dependence: the
Pearson $\rho$, the Spearman $\rho_{S}$ and the Kendall $\tau$ correlation
coefficients.

The rate of dependence between the variables $X$ and $Y$ is most often
measured by the Pearson correlation coefficient

$\rho=\frac{E(X-EX)(Y-EY)}{\sigma_{X}\sigma_{Y}}.$

The random variable

$\rho_{n}=\frac{\sum\limits_{i=1}^{n}(X_{i}-\bar{X})(Y_{i}-\bar{Y})}{\sqrt{\sum\limits_{i=1}^{n}(X_{i}-\bar{X})^{2}\sum\limits_{i=1}^{n}(Y_{i}-\bar{Y})^{2}}}$

is known as the sample Pearson correlation coefficient. It is a good
approximation for $\rho$ for large $n$ since $\rho_{n}\stackrel{{\scriptstyle
p}}{{\rightarrow}}\rho$; see, for example, Fisher (1921).

Let $X_{(1)}\leq\ldots\leq X_{(n)}$ be the order statistics obtained from the
sample $X_{1},\ldots,X_{n}$. For these order statistics, let us define their
concomitants $Y_{[1]},\ldots,Y_{[n]}$. Let $X_{i}=X_{(j)}$, then
$Y_{[j]}=Y_{i}$ be the concomitant of the order statistic $X_{(j)}$. The
concept of concomitants was proposed and discussed in the papers of David and
Galambos (1974), Bhattacharya (1974, 1984), Egorov and Nevzorov (1984), Goel
and Hall (1994), David (1994), Chu et al. (1999), David and Nagaraja (2003),
Balakrishnan and Lai (2009), Bairamov and Stepanov (2010), Balakrishnan and
Stepanov (2015). See also the references therein. Let

$R_{j}=\sum_{i=1}^{j-1}I(Y_{[i]}\leq Y_{[j]})\quad(1\leq i<j\leq n),$

where $I$ is the indicator function, be the rank of $Y_{[j]}$ amongst the
concomitants $Y_{[1]},\ldots,Y_{[j]}$. The random variable

$\tau_{n}=\frac{4\sum_{j=2}^{n}R_{j}}{n(n-1)}-1$

is known as the rank Kendall correlation coefficient. The theoretical analogue
of $\tau_{n}$ is presented by

$\tau=4E[F(X,Y)]-1=4\int_{\mathbb{R}^{2}}F(x,y)f(x,y)dxdy-1.$

It is known, see, for example, Balakrishnan and Lai (2009) or Stepanov (2024),
that $\tau=E\tau_{n}$ and $\tau_{n}\stackrel{{\scriptstyle
p}}{{\rightarrow}}\tau$. Let

$R_{j,n}=\sum_{i=1}^{n}I(Y_{[i]}\leq Y_{[j]})\quad(1\leq j\leq n)$

be the rank of $Y_{[j]}$ amongst the concomitants $Y_{[1]},\ldots,Y_{[n]}$.
The random variable

$\rho_{n,S}=1-\frac{6\sum_{j=1}^{n}(R_{j,n}-j)^{2}}{n^{3}-n}$

is known as the rank Spearman correlation coefficient. The theoretical
analogue of $\rho_{n,S}$ is presented by

$\rho_{S}=12E[H(X)G(Y)]-3=12\int_{\mathbb{R}^{2}}H(x)G(y)f(x,y)dxdy-3.$

The random variable $F(X,Y)-H(X)G(Y)$ is called the measure of quadrant
dependence; see, for example, p. 158 of Balakrishnan and Lai (2009).

The introduced above correlation coefficients $\rho,\ \rho_{S}$ and $\tau$ as
well as their statistical analogues $\rho_{n},\ \rho_{n,S}$ and $\tau_{n}$
have some advantages and disadvantages. Which correlation coefficient better
measures the dependence rate? For example, it is known that in many cases
$\mid\tau\mid\leq\mid\rho_{S}\mid\leq\mid\rho\mid$ and
$\mid\tau_{n}\mid\leq\mid\rho_{n,S}\mid\leq\mid\rho_{n}\mid$. Does it mean
that, say, $\rho$ is better than $\tau$? It is also known, for example, that
$\rho_{n}$ is very sensitive to contamination of outliers, while $\rho_{n,S}$
and $\tau_{n}$ are rather robust. This property can be viewed as an obvious
disadvantage of $\rho_{n}$ and, correspondingly, of $\rho$. Another
disadvantage of $\rho$ and, correspondingly of $\rho_{n}$, is that $\rho$
cannot be defined if the second moments do not exist. Observe that $\rho_{S}$
and $\tau$ exist for any absolutely continuous non-degenerate distribution
function. Another possible disadvantage of $\rho$ is that it is defined only
in terms of the first and second moments and, accordingly, loses much
information about the rate of dependence between $X$ and $Y$. On the other
hand, this can be viewed as an advantage since one can find $\rho$ knowing
only the first and second moments and does not need to know $F$. Since
$\tau_{n}$ and $\rho_{n,S}$ measure the dependence rate between ranks, rather
than between actual values of $X$ and $Y$, the correlation coefficients $\tau$
and $\rho_{S}$ are unaffected by any increasing transformation of $X$ and $Y$,
for example, by the monotone transformation based on the marginals of $F$.
Observe that $\rho$ is not invariant under this transformation. Then, if a
distribution function $F$ depends on the parameter $t$, and if this parameter
can be suppressed by one of the marginal transformations $u=H(x)$ or $v=G(y)$,
then $\rho_{S}$ and $\tau$ are free of $t$ for this distribution function $F$,
wherein $\rho$ can depend on $t$. In that case, the coefficients $\rho_{S}$
and $\tau$ are not elastic and, it seems that $\rho(t)$, which is elastic
here, measures the dependence rate better; see Example 3.3 below.

Although it is customary in bivariate data analysis to compute the dependence
rate by using only one correlation coefficient, in this paper we propose to
use a new correlation coefficient, which is a linear combination of $\rho_{S}$
and $\tau$. In general, before applying any correlation coefficient, a
researcher should study the scatter plot of data and decide which correlation
coefficient can be most helpful in that situation. The coefficient $\rho_{n}$
can be better used when the scatter plot indicates a linear relationship
between $X$ and $Y$, and $\rho_{n,S}$ and $\tau_{n}$ can be exploited in other
cases.

The rest of our paper is organized as follows. In Section 2, we propose a new
correlation coefficient, which is based on Kendal’s and Spearman’s correlation
coefficients. There, we also give a new extension of the Pearson correlation
coefficient. Section 3 contains examples, which illustrate the work of all
four correlation coefficients that we discuss in our paper. In the end of this
paper we present a conclusion.

##  2 New Correlation Coefficient. Extension of Pearson Correlation
Coefficient.

As was already pointed out, the following question is interesting: which
correlation coefficient measures the dependence rate better? Before answering
it, we analyze the known coefficients $\rho,\ \rho_{S}$ and $\tau$ and
introduce a new one $r$.

Let us consider the dependence function $F(x,y)-H(x)G(y)$, which gives us the
dependence value in every point $(x,y)$ of the distribution support. We adjust
the definition of the dependence function and redefine it as
$D(x,y)=6(F(x,y)-H(x)G(y))$. We then can view any correlation coefficient as
some average values of the dependence rate, or of the dependence function $D$.
Observe that the dependence function does not take into account the likelihood
of obtaining points $(x,y)$. Motivated by the above argument, we assume that
$r=E(D(X,Y))$ represents the “ideal” average of the dependence function.

New Correlation Coefficient In this paper, we propose to consider

$\displaystyle r$ $\displaystyle=$ $\displaystyle E(D(X,Y))$ $\displaystyle=$
$\displaystyle
6E[F(X,Y)-H(X)G(Y)]=6\int_{\mathbb{R}^{2}}(F(x,y)-H(x)G(y))f(x,y)dxdy$

as a new theoretical correlation coefficient. Observe that
$r=\frac{3\tau-\rho_{S}}{2}$. The following two inequalities

$\frac{3\tau-1}{2}\leq\rho_{S}\ (\mbox{if}\ \tau\geq
0),\quad\rho_{S}\leq\frac{1+3\tau}{2}\ (\mbox{if}\ \tau\leq 0)$

are given in (4.18) of Balakrishnan and Lai (2009). The first inequality
implies that $r\leq 1$, when the second one assumes that $r\geq-1$. Since in a
statistical experiment we usually do not know the underlying $F$, we need to
propose a statistical analogue of the correlation coefficient $r$. Taking into
account the definition of $r$, we propose to view

$r_{n}=1.5\tau_{n}-0.5\rho_{n,S}$

as the corresponding rank correlation coefficient. It is possible that $r_{n}$
can be presented in a more simple form. It was mentioned above that
$E\tau_{n}=\tau$. On page 157 of Balakrishnan and Lai (2009) one can find that
$E\rho_{n,S}=\frac{(n-2)\rho_{S}+3\tau}{n+1}$. It follows that

$Er_{n}=\frac{3\tau}{2}-\frac{(n-2)\rho_{S}+3\tau}{2(n+1)}=\frac{3n\tau-(n-2)\rho_{S}}{2(n+1)}\not=1.5\tau-0.5\rho_{S}=r.$

That way, although $r_{n}$ is a biased estimator of $r$, it is asymptotically
unbiased. The central limit theorem for $\rho_{n,S}$ was proved in Cifarelli
et al. (1996), from which we gain that $\rho_{n,S}\stackrel{{\scriptstyle
p}}{{\rightarrow}}\rho_{S}$. It is also known, see, for example, Stepanov
(2024), that $\tau_{n}\stackrel{{\scriptstyle p}}{{\rightarrow}}\tau\
(n\rightarrow\infty)$. By Slutsky’s theorem, we obtain that
$r_{n}\stackrel{{\scriptstyle p}}{{\rightarrow}}r$. We see that $r_{n}$ can be
considered as a statistical analogue of $r$.

Since $r$ is a mixture of $\rho_{S}$ and $\tau$, it inherits all advantages
and disadvantages of them. It is rather robust to the influence of outliers,
but it is constant when $F$ depends on a parameter $t$ and when $t$ can be
suppressed by the marginal transformations. Another situation when $r$
presumably works unsatisfactory (in comparison with $\rho$) is when the
dependence rate is close to $\pm 1$. By the definition of $r$, it is fair to
compare $r_{n}$ only with $\tau_{n}$ and $\rho_{n,S}$, because these three
coefficients measure dependence between the ranks. As was pointed out, $r$
represents some kind of the “ideal” average of the dependence function $D$. By
this argument, $r$ is better than $\rho_{S}$ and $\tau$. The correlation
coefficient $\rho$ represents the average of the random variable
$\frac{XY-E(XY)}{\sigma_{X}\sigma_{Y}}$. Accordingly, the coefficient $\rho$
is simpler, but some information about dependence can be lost, and $\rho$
cannot be so accurate when the situation is “far” from the linear relationship
between $X$ and $Y$. The coefficient $\rho_{n}$ measures dependence between
sample members and is somewhat different from the rank correlation
coefficients $\rho_{n,S},\ \tau_{n}$ and $r_{n}$, which measure association
between the ranks.

An interesting fact should be highlighted. By analyzing the dependence
function $D$, we conclude that the variables $X$ and $Y$ are independent iff
$Var(D(X,Y))=0$. We do not know how this discovery can be further used.

It should be noted that correlation coefficients can be compared by the sample
variances obtained by generation of these correlation coefficients. If a
sample variance is low, it indirectly indicates that the corresponding
correlation coefficient measures the dependence rate accurately, i.e., with
small squared error.

Let $\bar{\rho}_{n}$ and $S^{2}_{\rho_{n}}$ denote the mean and variance of a
sample consisting of generated Pearson correlation coefficients. Let the
designations for the other correlation coefficients are similar. Looking at
the results of simulation in Section 3, we find that generally

$|\bar{r}_{n}|\leq|\bar{\tau}_{n}|\leq|\bar{\rho}_{n,S}|\leq|\bar{\rho}_{n}|,$

and the inequalities

$S^{2}_{r_{n}}\leq S^{2}_{\tau_{n}}\leq S^{2}_{\rho_{n,S}}\leq
S^{2}_{\rho_{n}}$

hold true with the exception for the bivariate normal distribution, when $t$
is close to $\pm 1$. As we can see in Section 3 below, the rank correlation
coefficient $r_{n}$ inherits all the properties of $\rho_{S}$ and $\tau$, but
has one advantage over them. Most often it has the smallest sample variance.
Comparing the results of simulation in Section 3, we find that $\rho_{n,S},\
\tau_{n}$, and now $r_{n}$, “underestimate” the dependence rate in comparison
with $\rho_{n}$.

In the end of this section, we would like to propose a new extension of the
Pearson correlation coefficient $\rho$. In the case when $Var(X)$ (or
$Var(Y)$) does not exist, the coefficient $\rho$ is formally not defined. It
was known that it can be generalized and redefined, for example, for the
bivariate stable distribution functions; see p. 146 of Balakrishnan and Lai
(2007). Thus, for the symmetric bivariate Cauchy distribution function, we can
apply the Cauchy principal value. Then, for this distribution function,
$cov(X,Y)=0$ and we can put $\rho=0$.

In this work, we propose a new extension of $\rho$. In the case, when the
second moments do not exist, but simulation experiments show that
$\bar{\rho}_{n}\approx\tilde{\rho}$ for large $n$, we propose to consider
$\tilde{\rho}$ as a generalization of $\rho$ provided that the values of
$\bar{\rho}_{n}$ are “close” to the values of $\bar{\rho}_{n,S},\
\bar{\tau}_{n}$ and $\bar{r}_{n}$; see Examples 3.2 below.

##  3 Examples

Example 3.1 Let

$F_{t}(x,y)=\frac{1}{2\pi\sqrt{1-t^{2}}}\int_{-\infty}^{x}\int_{-\infty}^{y}e^{-\frac{u^{2}-2tuv+v^{2}}{2(1-t^{2})}}dudv\quad(x,y\in\mathbb{R},\
-1<t<1)$

be a bivariate normal distribution. We have

$\rho(t)=t,\quad\rho_{S}(t)=\frac{6}{\pi}asin(t/2),\quad\tau(t)=\frac{2}{\pi}asin(t)\quad
r(t)=\frac{3}{\pi}\left(asin(t)-asin(t/2)\right).$

For visual comparison, we present the above correlation coefficients for
different values of $t$ in Table 3.1.

Table 3.1 Values of correlation coefficients for different $t$.

| $t=-0.99$ | $t=-0.70$ | $t=-0.30$ | $t=-0.10$ | $t=0.10$ | $t=0.30$ | $t=0.70$ | $t=0.99$
---|---|---|---|---|---|---|---|---
$\rho$ | -0.9900 | -0.7000 | -0.3000 | -0.1000 | 0.1000 | 0.3000 | 0.7000 | 0.9900
$\rho_{s}$ | -0.9890 | -0.6829 | -0.2876 | -0.0955 | 0.0955 | 0.2876 | 0.6829 | 0.9890
$\tau$ | -0.9098 | -0.4936 | -0.1940 | -0.0638 | 0.0638 | 0.1940 | 0.4936 | 0.9098
$r$ | -0.8703 | -0.3990 | -0.1472 | -0.0479 | 0.0479 | 0.1472 | 0.3990 | 0.8703

It is clear from Table 3.1 that $\mid
r|\leq|\tau|\leq|\rho_{S}|\leq|\rho\mid$. We conducted several simulation
experiments. We generated one thousand times these correlation coefficients
for samples $X_{1},\ldots,X_{1000}$ and $Y_{1},\ldots,Y_{1000}$, computed
their means and variances and presented them in Tables 3.2–3.3, respectively.

Table 3.2 Sample means of correlation coefficients for different $t$.

| $t=-0.99$ | $t=-0.70$ | $t=-0.30$ | $t=-0.10$ | $t=0.10$ | $t=0.30$ | $t=0.70$ | $t=0.99$
---|---|---|---|---|---|---|---|---
$\bar{\rho}_{n}$ | -0.9900 | -0.7000 | -0.3002 | -0.0999 | 0.0996 | 0.3003 | 0.6985 | 0.9900
$\bar{\rho}_{n,S}$ | -0.9887 | -0.6825 | -0.2874 | -0.0955 | 0.0949 | 0.2878 | 0.6812 | 0.9888
$\bar{\tau}_{n}$ | -0.9100 | -0.4940 | -0.1940 | -0.0638 | 0.0634 | 0.1944 | 0.4926 | 0.9100
$\bar{r}_{n}$ | -0.8706 | -0.3994 | -0.1473 | -0.0480 | 0.0477 | 0.1476 | 0.3982 | 0.8706

Table 3.3 Sample variances of correlation coefficients for different $t$.

| $t=-0.99$ | $t=-0.70$ | $t=-0.30$ | $t=-0.10$ | $t=0.10$ | $t=0.30$ | $t=0.70$ | $t=0.99$
---|---|---|---|---|---|---|---|---
$S^{2}_{\rho_{n}}$ | $4.1\cdot 10^{-7}$ | $2.9\cdot 10^{-4}$ | $8.1\cdot 10^{-4}$ | $1.0\cdot 10^{-3}$ | $9.0\cdot 10^{-4}$ | $8.4\cdot 10^{-4}$ | $2.6\cdot 10^{-4}$ | $4.0\cdot 10^{-7}$
$S^{2}_{\rho_{n,S}}$ | $7.0\cdot 10^{-7}$ | $3.6\cdot 10^{-4}$ | $8.5\cdot 10^{-4}$ | $1.0\cdot 10^{-3}$ | $9.0\cdot 10^{-4}$ | $8.8\cdot 10^{-4}$ | $3.2\cdot 10^{-4}$ | $7.2\cdot 10^{-7}$
$S^{2}_{\tau_{n}}$ | $9.9\cdot 10^{-6}$ | $2.6\cdot 10^{-4}$ | $4.0\cdot 10^{-4}$ | $4.5\cdot 10^{-4}$ | $4.1\cdot 10^{-4}$ | $4.2\cdot 10^{-4}$ | $2.3\cdot 10^{-4}$ | $9.8\cdot 10^{-6}$
$S^{2}_{r_{n}}$ | $1.8\cdot 10^{-5}$ | $2.2\cdot 10^{-4}$ | $2.4\cdot 10^{-4}$ | $2.5\cdot 10^{-4}$ | $2.3\cdot 10^{-4}$ | $2.6\cdot 10^{-4}$ | $1.9\cdot 10^{-4}$ | $1.8\cdot 10^{-5}$

It seems that $r_{n}$ is a good rank correlation coefficient, since it has the
uniformly small sample variance.

We generated once more one thousand times correlation coefficients $\rho_{n},\
\rho_{n,S},\ \tau_{n}$ and $r_{n}$. One thousand observations were taken from
$F$ and five observations were outliers. We obtained the following results for
$t=-0.99$

$\bar{\rho}_{n}=-0.3305,\ \bar{\rho}_{n,S}=-0.9592,\ \bar{\tau}_{n}=-0.8910,\
\bar{r}_{n}=-0.8569.$

As we can see, the effect of outliers on $\rho_{n}$ is great.

We also simulated the correlation coefficients that measure correlation
between the correlation coefficients themselves. Let, for example,
$\rho_{100}^{(i)},\rho_{100,S}^{(j)}\ (1\leq i,j\leq 100)$ be samples of size
100 of Pearson’s and Spearman’s correlation coefficients and
$\bar{\rho}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})$ be the mean of
Pearson’s correlation coefficients, which measure the dependence rate between
the samples consisting of Pearson’s and Spearman’s correlation coefficients.
For $t=-0.99$, we present all results for all types of correlation
coefficients:

$\bar{\rho}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=0.7340,\quad\bar{\rho}_{100}(\rho_{100}^{(i)},\tau_{100}^{(j)})=0.7987,\quad\bar{\rho}_{100}(\rho_{100}^{(i)},r_{100}^{(j)})=0.8017,$
$\bar{\rho}_{100}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9578,\quad\bar{\rho}_{100}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.9485,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9985,$
$\bar{\rho}_{100,S}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=0.7242,\quad\bar{\rho}_{100,S}(\rho_{100}^{(i)},\tau_{100}^{(j)})=0.7826,\quad\bar{\rho}_{100,S}(\rho_{100}^{(i)},r_{100}^{(j)})=0.7856,$
$\bar{\rho}_{100,S}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9568,\quad\bar{\rho}_{100,S}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.9479,\quad\bar{\rho}_{100,S}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9992,$
$\bar{\tau}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=0.5385,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},\tau_{100}^{(j)})=0.5972,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100}^{(j)})=0.5973,$
$\bar{\tau}_{100}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9578,\quad\bar{\tau}_{100}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.8121,\quad\bar{\tau}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9861.$

As we can see the correlations between the rank correlation coefficients are
rather close to one. The Pearson sample correlation coefficient is not so much
correlated with the rank correlation coefficients. The smallest and largest
correlation coefficients for the normal case for $t=-0.99$ are
$\bar{\tau}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=0.5385$ and
$\bar{\rho}_{100,S}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9992$, respectively. The
smallest and largest correlation coefficients for $t=-0.7$ are,
correspondingly,
$\bar{\tau}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=0.7126$ and
$\bar{\rho}_{100,S}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9956$. As $t$ advances
to zero the difference between the largest and smallest correlation
coefficients decreases. Thus, for $t=-0.1$, we have
$\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100}^{(j)})=0.7786$ and
$\bar{\rho}_{100,S}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9956$. For positive
$t$ the situation is symmetric. We conclude that there is some similarity in
behaviors of $\rho_{n}$ and the rank correlation coefficients when $t$ is not
close to $\pm 1$.

Example 3.2 Let us consider the bivariate Pareto distribution function

$F_{t}(x,y)=1-\frac{1}{(y+1)^{t}}-\frac{1}{(x+1)^{t}}+\frac{1}{(x+y+1)^{t}}\quad(x>0,\
y>0,\ t>0).$

We have $\rho(t)=\frac{1}{t}\ (t>2)$ and $\tau(t)=\frac{1}{2t+1}\ (t>0)$.
Since correlation coefficients $\rho_{s}(t)$ and $r(t)$ can be found only
numerically, we present in Table 3.4 their values for different $t$ along with
the values of $\rho(t)$ and $\tau(t)$. The means and variances of these
correlation coefficients are presented in Tables 3.5–3.6, respectively.

Table 3.4 Values of correlation coefficients for different $t$.

| $t=0.1$ | $t=0.5$ | $t=1$ | $t=2.1$ | $t=5$ | $t=10$ | $t=50$ | $t=100$
---|---|---|---|---|---|---|---|---
$\rho$ | – | – | – | 0.4761 | 0.2000 | 0.1000 | 0.0200 | 0.0100
$\rho_{s}$ | 0.92368 | 0.6822 | 0.4784 | 0.2839 | 0.1358 | 0.0714 | 0.0149 | 0.0075
$\tau$ | 0.8333 | 0.5000 | 0.3333 | 0.1923 | 0.0909 | 0.0476 | 0.0099 | 0.0050
$r$ | 0.6903 | 0.4089 | 0.2608 | 0.1465 | 0.0684 | 0.0357 | 0.0074 | 0.0037

Table 3.5 Sample means of correlation coefficients for different $t$.

| $t=0.1$ | $t=0.5$ | $t=1$ | $t=2.1$ | $t=5$ | $t=10$ | $t=50$ | $t=100$
---|---|---|---|---|---|---|---|---
$\bar{\rho}_{n}$ | 0.9545 | 0.7972 | 0.6502 | 0.4401 | 0.1996 | 0.1003 | 0.0198 | 0.0099
$\bar{\rho}_{n,S}$ | 0.9580 | 0.6818 | 0.4784 | 0.2831 | 0.1364 | 0.0709 | 0.0143 | 0.0073
$\bar{\tau}_{n}$ | 0.8334 | 0.5000 | 0.3337 | 0.1919 | 0.0914 | 0.0474 | 0.0096 | 0.0049
$\bar{r}_{n}$ | 0.7712 | 0.4091 | 0.2613 | 0.1463 | 0.0689 | 0.0356 | 0.0072 | 0.0037

Table 3.6 Sample variances of correlation coefficients for different $t$.

| $t=0.1$ | $t=0.5$ | $t=1$ | $t=2.1$ | $t=5$ | $t=10$ | $t=50$ | $t=100$
---|---|---|---|---|---|---|---|---
$S^{2}_{\rho_{n}}$ | $1.9\cdot 10^{-2}$ | $5.4\cdot 10^{-2}$ | $5.2\cdot 10^{-2}$ | $2.2\cdot 10^{-2}$ | $3.5\cdot 10^{-3}$ | $1.6\cdot 10^{-3}$ | $9.5\cdot 10^{-4}$ | $1.0\cdot 10^{-3}$
$S^{2}_{\rho_{n,S}}$ | $1.4\cdot 10^{-5}$ | $4.0\cdot 10^{-4}$ | $7.1\cdot 10^{-4}$ | $9.6\cdot 10^{-4}$ | $11.0\cdot 10^{-4}$ | $9.7\cdot 10^{-4}$ | $9.9\cdot 10^{-4}$ | $1.0\cdot 10^{-3}$
$S^{2}_{\tau_{n}}$ | $4.9\cdot 10^{-5}$ | $2.9\cdot 10^{-4}$ | $3.9\cdot 10^{-4}$ | $4.7\cdot 10^{-4}$ | $5.0\cdot 10^{-4}$ | $4.3\cdot 10^{-4}$ | $4.4\cdot 10^{-4}$ | $4.5\cdot 10^{-4}$
$S^{2}_{r_{n}}$ | $7.6\cdot 10^{-5}$ | $2.5\cdot 10^{-4}$ | $2.7\cdot 10^{-4}$ | $2.8\cdot 10^{-4}$ | $2.9\cdot 10^{-4}$ | $2.5\cdot 10^{-4}$ | $2.5\cdot 10^{-4}$ | $2.5\cdot 10^{-4}$

It seems that again $r$ has the smallest sample variance and can be viewed as
a solid correlation coefficient.

It should be noted that $\rho_{n}$ measures the dependence rate when
$t\in(0,2]$ and when its theoretical analogue $\rho$ does not exist. Moreover,
it follows the trend given by the rank correlation coefficients, i.e.,
$\bar{\rho}_{n}(t)$ tends to one as $t$ tends to zero. Taking into account the
last paragraph of introduction, we define the correlation coefficient
$\tilde{\rho}$, which is an extension of $\rho$. Let $\tilde{\rho}$ be such
that $\tilde{\rho}(t)=\rho(t)$ if $t>2$ and $\tilde{\rho}(t)=\bar{\rho}_{n}$,
where $\bar{\rho}_{n}(t)$ is found by simulation for every $t\in(0,2]$.

We also simulated the correlation coefficients which measure correlation
between the correlation coefficients themselves. For $t=0.1$ we present below
all results for all types of correlation coefficients:

$\bar{\rho}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=-2.8\cdot
10^{-3},\quad\bar{\rho}_{100}(\rho_{100}^{(i)},\tau_{100}^{(j)})=-3.9\cdot
10^{-4},\quad\bar{\rho}_{100}(\rho_{100}^{(i)},r_{100}^{(j)})=1.5\cdot
10^{-4},$
$\bar{\rho}_{100}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9655,\quad\bar{\rho}_{100}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.9489,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9983,$
$\bar{\rho}_{100,S}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=-4.2\cdot
10^{-3},\quad\bar{\rho}_{100,S}(\rho_{100}^{(i)},\tau_{100}^{(j)})=-2.4\cdot
10^{-3},\quad\bar{\rho}_{100,S}(\rho_{100}^{(i)},r_{100}^{(j)})=-2.8\cdot
10^{-3},$
$\bar{\rho}_{100,S}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9646,\quad\bar{\rho}_{100,S}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.9491,\quad\bar{\rho}_{100,S}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9976,$
$\bar{\tau}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=-2.7\cdot
10^{-3},\quad\bar{\tau}_{100}(\rho_{100}^{(i)},\tau_{100}^{(j)})=-1.4\cdot
10^{-3},\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100}^{(j)})=-1.6\cdot
10^{-3},$
$\bar{\tau}_{100}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9655,\quad\bar{\tau}_{100}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.8132,\quad\bar{\tau}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9693.$

As we can see again, the correlations between the rank correlation
coefficients are rather close to one. On the contrary, $\rho_{n}$ is almost
uncorrelated with the rank correlation coefficients, though the value of
$\bar{\rho}_{n}$ is rather close to the values of $\bar{\rho}_{n,S},\
\bar{\tau}_{n}$ and $\bar{r}_{n}$. When $t$ increases the situation becomes
more even. For $t=1,10,100,1000$ we present only the smallest and largest
correlation coefficients:

$t=1,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=0.0358,\quad\bar{\tau}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9961,$
$t=10,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100,S}^{(j)})=0.4558,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9977,$
$t=100,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100,S}^{(j)})=0.5507,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9977.$
$t=1000,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100,S}^{(j)})=0.5607,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9979.$

As $t$ increases further, the difference between the largest and smallest
correlation coefficients does not change much.

Example 3.3 Let

$F_{t}(x,y)=1-e^{-x}-\frac{1}{(1+y)^{t}}+\frac{e^{-x(1+y)^{t}}}{(1+y)^{t}}\quad(x,y>0).$

We have $\rho=-\frac{\sqrt{t(t-2)}}{2t-1}\ (t>2)$. For $t>0$ the rank
correlation coefficients are free of $t$ and take the following values

$\rho_{S}=12\log 2-9=-0.6822,\quad\tau=-0.5,\quad r=15/4-6\log 2=-0.4089.$

One can see that here only $\rho$ depends on $t$ and presumably better
reflects the dependence rate when $t>2$. The values of $\rho$ along with the
means of all correlation coefficients are presented in Table 3.7 for $t>0$.
The corresponding sample variances are presented in Table 3.8.

Table 3.7 Values of $\rho$ and sample means of correlation coefficients.

| $t=0.1$ | $t=0.5$ | $t=1$ | $t=2.1$ | $t=5$ | $t=10$ | $t=50$ | $t=100$
---|---|---|---|---|---|---|---|---
$\rho$ | – | – | – | -0.1432 | -0.4303 | -0.4708 | -0.4948 | -0.4975
$\bar{\rho}_{n}$ | -0.0324 | -0.0525 | -0.1192 | -0.2998 | -0.4364 | -0.4726 | -0.4961 | -0.4989
$\bar{\rho}_{n,S}$ | -0.6820 | -0.6816 | -0.6814 | -0.6819 | -0.6815 | -0.6811 | -0.6814 | -0.6820
$\bar{\tau}_{n}$ | -0.5002 | -0.4999 | -0.4998 | -5.002 | -0.4999 | -0.4995 | -0.4998 | -0.5002
$\bar{r}_{n}$ | -0.4093 | -0.4091 | -0.4090 | -0.4093 | -0.4091 | -0.4087 | -0.4089 | -0.4093

Table 3.8 Sample variances of correlation coefficients for different $t$.

| $t=0.1$ | $t=0.5$ | $t=1$ | $t=2.1$ | $t=5$ | $t=10$ | $t=50$ | $t=100$
---|---|---|---|---|---|---|---|---
$S^{2}_{\rho_{n}}$ | $2.5\cdot 10^{-5}$ | $3.7\cdot 10^{-4}$ | $2.7\cdot 10^{-3}$ | $4.2\cdot 10^{-3}$ | $7.8\cdot 10^{-4}$ | $3.7\cdot 10^{-4}$ | $2.5\cdot 10^{-4}$ | $2.5\cdot 10^{-4}$
$S^{2}_{\rho_{n,S}}$ | $3.7\cdot 10^{-4}$ | $3.7\cdot 10^{-4}$ | $3.8\cdot 10^{-4}$ | $3.8\cdot 10^{-4}$ | $3.8\cdot 10^{-4}$ | $3.8\cdot 10^{-4}$ | $3.7\cdot 10^{-4}$ | $3.8\cdot 10^{-4}$
$S^{2}_{\tau_{n}}$ | $2.7\cdot 10^{-4}$ | $2.7\cdot 10^{-4}$ | $2.7\cdot 10^{-4}$ | $2.7\cdot 10^{-4}$ | $2.7\cdot 10^{-4}$ | $2.7\cdot 10^{-4}$ | $2.6\cdot 10^{-4}$ | $2.7\cdot 10^{-4}$
$S^{2}_{r_{n}}$ | $2.3\cdot 10^{-4}$ | $2.3\cdot 10^{-4}$ | $2.3\cdot 10^{-4}$ | $2.3\cdot 10^{-4}$ | $2.3\cdot 10^{-4}$ | $2.3\cdot 10^{-4}$ | $2.2\cdot 10^{-4}$ | $2.3\cdot 10^{-4}$

The coefficient $\rho_{n}$ again measures the dependence rate for $t\in(0,2]$
when its theoretical analogue $\rho$ is not defined. However here, $\rho_{n}$
does not follow the trend of the rank correlation coefficients. Observe that
for $t=0$, we have

$P(X\leq x,Y\leq y)=P(X\leq x)\cdot 1\quad\mbox{and}\quad P(Y\leq y)=0.$

That way, the case $t=0$ is not associated with independence, but $\rho_{n}$
wrongly indicates it. In general, the coefficient $\rho_{n}$ gives us wrong
information when $t\in(0,2]$ in spite of the statistic $S^{2}_{\rho_{n}}$
takes there on small values. The values of $\bar{\rho}_{n}$ are small for
$t\in(0,2]$, because for $t=0.1$ we have $cov_{sam}(X,Y)=-7.6*10\cdot^{30},\
S_{x}=0.9542$ and $S_{y}=2.5\cdot 10^{32}$, where $cov_{sam}$ is the sample
covariance. The independence implies that the data is wildly dispersed in both
directions $x$ and $y$. In this example it happens only in one direction,
along the variable $y$.

We also simulated the correlation coefficients between correlation
coefficients themselves. For $t=0.1$ we present below all results for all
types of correlation coefficients:

$\bar{\rho}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=-0.1265,\quad\bar{\rho}_{100}(\rho_{100}^{(i)},\tau_{100}^{(j)})=-0.1334,\quad\bar{\rho}_{100}(\rho_{100}^{(i)},r_{100}^{(j)})=-0.1360,$
$\bar{\rho}_{100}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9873,\quad\bar{\rho}_{100}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.9660,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9948,$
$\bar{\rho}_{100,S}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=-0.1504,\quad\bar{\rho}_{100,S}(\rho_{100}^{(i)},\tau_{100}^{(j)})=-0.1574,\quad\bar{\rho}_{100,S}(\rho_{100}^{(i)},r_{100}^{(j)})=-0.1601,$
$\bar{\rho}_{100,S}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9842,\quad\bar{\rho}_{100,S}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.9609,\quad\bar{\rho}_{100,S}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9929,$
$\bar{\tau}_{100}(\rho_{100}^{(i)},\rho_{100,S}^{(j)})=-0.1021,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},\tau_{100}^{(j)})=-0.1068,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100}^{(j)})=-0.1088,$
$\bar{\tau}_{100}(\rho_{100,S}^{(i)},\tau_{100}^{(j)})=0.9873,\quad\bar{\tau}_{100}(\rho_{100,S}^{(i)},r_{100}^{(j)})=0.8393,\quad\bar{\tau}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9382.$

As we can see, the correlations between the rank correlation coefficients are
rather close to one. The coefficient $\rho_{n}$ is even negatively correlated
with the rank correlation coefficients.

When $t$ increases the situation does not change much. For
$t=1,10,100,1000,10000$, we present the smallest and largest correlation
coefficients:

$t=1,\quad\bar{\rho}_{100,S}(\rho_{100}^{(i)},\rho_{100}^{(j)})=-0.0514,\quad\bar{\tau}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9947.$
$t=2,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100}^{(j)})=0.0634,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9948.$
$t=10,\quad\bar{\rho}_{100,S}(\rho_{100}^{(i)},r_{100,S}^{(j)})=0.4376,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9948.$
$t=100,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100,S}^{(j)})=0.3870,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9948.$
$t=1000,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100,S}^{(j)})=0.4037,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9948.$
$t=10000,\quad\bar{\tau}_{100}(\rho_{100}^{(i)},r_{100,S}^{(j)})=0.3903,\quad\bar{\rho}_{100}(\tau_{100}^{(i)},r_{100}^{(j)})=0.9948.$

As was pointed out in Introduction, the coefficients $\tau$ and $\rho_{S}$,
and then, consequently, $r$ are invariant under any increasing transformation.
Let $1-\frac{1}{(1+y)^{t}}=v\in[0,1]$. Then

$F_{1}(x,v)=v-e^{-x}+(1-v)e^{-\frac{x}{1-v}}\quad(x\geq 0,\ 0\leq v\leq 1).$

For $F_{1}$, we have the same results for $\rho_{S},\ \tau$ and $r$ as we had
earlier for $F$, i.e.,

$\rho=-\frac{1}{\sqrt{3}}=-0.5774,\quad\rho_{S}=12\log
2-9=-0.6822,\quad\tau=-0.5,\quad r=15/4-6\log 2=-0.4089.$

We present the means and variances of one thousand simulations of these
correlation coefficients computed for samples $X_{1},\ldots,X_{1000}$ and
$Y_{1},\ldots,Y_{1000}$ with $F_{1}$:

$\bar{\rho}_{n}=-0.5780,\quad\bar{\rho}_{n,S}=-0.6816,\quad\bar{\tau}_{n}=-0.4999,\quad\bar{r}_{n}=-0.4091,$
$S^{2}_{\rho_{n}}=3.2\cdot 10^{-4},\quad S^{2}_{\rho_{n,S}}=3.8\cdot
10^{-4},\quad S^{2}_{\tau_{n}}=2.7\cdot 10^{-4},\quad S^{2}_{r_{n}}=2.3\cdot
10^{-4}.$

Example 3.4 Let

$F(x,y)=\frac{x(y+1)(1+\alpha(1-x)y)}{2}\quad(0\leq x\leq 1,\ -1\leq y\leq 0)$
$F(x,y)=\frac{x}{2}+\frac{xy(1+\alpha(1-x)(1-y))}{2}\quad(0\leq x\leq 1,\
0<y\leq 1),$

where $t\in[-1,1]$. If $t\not=0$, the variables $X$ and $Y$ are dependent, but

$\tau=\tau_{S}=\rho=r=0.$

Between uncorrelatedness and independence lies semi-independence that is
defined by the two equalities $E(Y|X)=EY$ and $E(X|Y)=EX$; see Jensen (1988).
However here, $E(Y|X)=EY=0$, but

$\frac{1}{2}-\frac{\alpha}{6}=E(X|Y)\not=EX=\frac{1}{2}.$

The corresponding simulation results are presented in Tables 3.9–3.10 below.

Table 3.9 Sample means of correlation coefficients for different $t$.

| $t=-1$ | $t=-0.5$ | $t=-0.1$ | $t=0$ | $t=0.1$ | $t=0.5$ | $t=1$
---|---|---|---|---|---|---|---
$\bar{\rho}_{n}$ | $-1.0\cdot 10^{-3}$ | $-2.1\cdot 10^{-4}$ | $8.6\cdot 10^{-4}$ | $-2.0\cdot 10^{-3}$ | $3.5\cdot 10^{-4}$ | $-8.5\cdot 10^{-4}$ | $1.6\cdot 10^{-3}$
$\bar{\rho}_{n,S}$ | $-1.3\cdot 10^{-3}$ | $-2.1\cdot 10^{-3}$ | $9.2\cdot 10^{-4}$ | $-2.0\cdot 10^{-3}$ | $3.0\cdot 10^{-4}$ | $-7.6\cdot 10^{-4}$ | $1.4\cdot 10^{-3}$
$\bar{\tau}_{n}$ | $-9.2\cdot 10^{-4}$ | $-1.4\cdot 10^{-3}$ | $6.1\cdot 10^{-4}$ | $-1.8\cdot 10^{-3}$ | $2.2\cdot 10^{-4}$ | $-4.9\cdot 10^{-4}$ | $9.9\cdot 10^{-4}$
$\bar{r}_{n}$ | $-7.2\cdot 10^{-4}$ | $-1.0\cdot 10^{-3}$ | $4.6\cdot 10^{-4}$ | $-1.0\cdot 10^{-3}$ | $1.7\cdot 10^{-4}$ | $-3.6\cdot 10^{-4}$ | $7.7\cdot 10^{-4}$

Table 3.10 Sample variances of correlation coefficients for different $t$.

| $t=-1$ | $t=-0.5$ | $t=-0.1$ | $t=0$ | $t=0.1$ | $t=0.5$ | $t=1$
---|---|---|---|---|---|---|---
$S^{2}_{\rho_{n}}$ | $9.7\cdot 10^{-4}$ | $9.4\cdot 10^{-4}$ | $9.9\cdot 10^{-4}$ | $9.4\cdot 10^{-4}$ | $9.7\cdot 10^{-4}$ | $1.1\cdot 10^{-3}$ | $9.9\cdot 10^{-4}$
$S^{2}_{\rho_{n,S}}$ | $1.0\cdot 10^{-3}$ | $9.5\cdot 10^{-4}$ | $9.9\cdot 10^{-4}$ | $9.4\cdot 10^{-4}$ | $9.7\cdot 10^{-4}$ | $1.1\cdot 10^{-3}$ | $1.0\cdot 10^{-3}$
$S^{2}_{\tau_{n}}$ | $4.6\cdot 10^{-4}$ | $4.2\cdot 10^{-4}$ | $4.4\cdot 10^{-4}$ | $4.2\cdot 10^{-4}$ | $4.3\cdot 10^{-4}$ | $4.7\cdot 10^{-4}$ | $4.7\cdot 10^{-4}$
$S^{2}_{r_{n}}$ | $2.6\cdot 10^{-4}$ | $2.4\cdot 10^{-4}$ | $2.5\cdot 10^{-4}$ | $2.4\cdot 10^{-4}$ | $2.4\cdot 10^{-4}$ | $2.6\cdot 10^{-4}$ | $2.7\cdot 10^{-4}$

The coefficient $r_{n}$ has the smallest sample variance and can again be
viewed as a good rank correlation coefficient. It is interesting that the
generation of correlation coefficients for $t=0.5$ gave us negative numbers.

## Conclusion

In the present paper, we have discussed the Pearson $\rho$, Spearman
$\rho_{S}$, Kendall $\tau$ correlation coefficients and have introduced a new
correlation coefficient $r$. We have proposed to consider $\tilde{\rho}$, such
that $\rho_{n}$ approximates $\tilde{\rho}$ in simulation experiments for
large $n$, as a new extension of $\rho$.

We have conducted simulation experiments and studied the behavior of the above
correlation coefficients. We have come to the conclusion that the behavior of
$\rho_{n}$ can be very different from the behavior of the rank correlation
coefficients $\rho_{n,S},\ \tau_{n},\ r_{n}$. This follows from the definition
of $\rho_{n}$ and $\rho_{n,S},\ \tau_{n}\ r_{n}$. The coefficient $\rho_{n}$
measures the association between variables $X_{i}$ and $Y_{j}$, whereas the
coefficients $\rho_{n,S},\ \tau_{n}$ and $r_{n}$ measure the association
between their ranks. In this paper we have raised the question: which
correlation coefficient measures the dependence rate better? We think that it
is fair to compare only $\tau_{n},\ \rho_{n,S}$ and $r_{n}$, since these three
coefficients are the ranks correlation coefficients. We have come to the
conclusion that $r$ looks a little better than $\rho_{S}$ and $\tau$. The
coefficient $\rho$, from one hand, and the coefficients $\rho_{S},\ \tau,\ r$,
from the other hand, measure in their own peculiar ways, and cannot be
compared generally. We recommend to use $\rho$ when the association between
$X$ and $Y$ is close to $\pm 1$ and to use $\rho_{S},\ \tau,\ r$, otherwise.

## References

Bairamov, I., Stepanov, A. (2010). Numbers of near-maxima for the bivariate
case, Statistics $\&$ Probability Letters, 80, 196–205.

Balakrishnan, N., Lai, C.D. (2009). Continuous Bivariate Distributions, Second
edition, Springer.

Balakrishnan, N., Stepanov, A. (2015). Limit results for concomitants of order
statistics, Metrika, 78, 385–397.

Bhattacharya, B.B. (1974). Convergence of sample paths of normalized sums of
induced order statistics, Ann. Statist. , 2, 1034–1039.

Bhattacharya, B.B. (1984). Induced order statistics: theory and applications,
In Handbook of Statistics 4, Ed. Krishnaiah, P. R., Sen, P. K., North Holland,
Amsterdam, 383–403.

Chu, S.J., Huang, W.J., Chen, H. (1999). A study of asymptotic distributions
of concomitants of certain order statistics, Statist. Sinica, 9, 811–830.

Cifarelli, D. M., Conti, P. L., Regazzini, E. (1996). On the asymptotic
distribution of a general measure of monotone dependence, The Annals of
Statistics, 25, 1386–1399.

Daniels, H. E. (1950). Rank correlation and population models, Journal of the
Royal Statistical Society, Ser. B, 12 (2), 171–191.

David, H.A. (1994). Concomitants of Extreme Order Statistics, In Extreme Value
Theory and Applications, Proceedings of the Conference on Extreme Value Theory
and Applications, 1, Ed. Galambos, J., Lechner, J., Simiu, E., Kluwer Academic
Publishers, Boston 211–224.

David, H.A., Galambos, J. (1974). The asymptotic theory of concomitants of
order statistics, J. Appl. Probab., 11, 762–770.

David, H.A., Nagaraja, H.N. (2003). Order Statistics, Third edition, John
Wiley & Sons, NY.

Durbin J., Stuart, A. (1951). Inversions and rank correlation coefficients,
Journal of the Royal Statistical Society, Ser. B 13 (2), 303–309.

Egorov, V. A., Nevzorov, V. B. (1984). Rate of convergence to the Normal law
of sums of induced order statistics, Journal of Soviet Mathematics (New York),
25, 1139–1146.

Goel, P. K., Hall, P. (1994). On the average difference between concomitants
and order statistics, Ann. Probab., 22, 126–144.

Fisher, R. A. (1921). On the ”probable error” of a coefficient of correlation
deduced from a small sample, Metron, 1, 3–32.

Jensen, D.R. (1988). Semi-independence. In: Encyclopedia of Statistical
Sciences, Volume 8, S. Kotz and N.L. Johnson (eds.), 358 -359. John Wiley and
Sons, New York?

Kendall, M. G. (1970). Rank Correlation Methods, London, Griffin.

Omey, E., Gulk, S.Van (2008). Central limit theorems for variances and
correlation coefficients, Preprint, doi: 10.13140/RG.2.2.21663.66727.

Stepanov, A.V. (2024). On Kendall’s correlation coefficient, Vestnik St.
Petersburg University, Mathematics, accepted.

Xu, W., Hou, Y., Hung, Y. S., Zou, Y. (2009). Comparison of Spearman’s rho and
Kendall’s tau in Normal and Contaminated Normal Models, arXiv:1011.2009v1
[cs.IT].


# Propagation gap for shear waves in binary liquids: Analytical and simulation
study

Taras Bryk1,2, Maria Kopcha1 1 Institute for Condensed Matter Physics of NAS
of Ukraine,
UA-79011 Lviv, Ukraine  2Institute of Applied Mathematics and Fundamental
Sciences,
Lviv National Polytechnic University, UA-79013 Lviv, Ukraine

###### Abstract

Transverse collective excitations in a binary liquid mixture are studied for
different mass ratio of components $R$ at fixed numerical density. Increasing
the mass ratio results in a larger gap width for shear waves. A four-variable
dynamic model of transverse dynamics in binary liquids is solved analytically
with an account of cross-correlations between total mass and mass-
concentration current fluctuations. An equation for the propagation gap of
shear waves for binary liquids is reported and analyzed.

##  I Introduction

Transverse dynamics in liquids is very simple on macroscopic scales, where the
only shear relaxation perfectly describes the long-time behavior of wave-
number-dependent transverse mass-current autocorrelations for sufficiently
small wave numbers $k$Boo ; Han :

$F^{T}_{JJ}(k,t)\stackrel{{\scriptstyle k\to 0}}{{=}}mk_{B}Te^{-D_{v}k^{2}t}$
(1)

with $D_{v}=\eta/\rho$ being kinematic viscosity, $\eta$ \- shear viscosity,
and $\rho$ \- mass density. This fact is a consequence of existing the only
transverse conserved quantity - transverse component of total momentum - that
excludes any transverse wave-like motion on macroscopic scale in liquid, i.e.
explains the absence of transverse sound in liquid systems. From (1) it
follows, that the hydrodynamic long-time behavior of shear relaxation for any
liquid (simple or complex) will be qualitatively similar. Outside the
hydrodynamic regime, with increasing wave numbers, the overdamped propagating
transverse excitations called shear waves can emerge in the liquid. Shear
waves can propagate only on the nanoscales because of the strong damping and
this fact makes them absolutely different from the transverse sound in
crystalline and glassy solids. Inability of the shear waves to exist on
macroscopic distances is reflected in their dispersion law $\omega(k)$, which
was obtained in the end of 1960-sChu69 ; Boo , and makes evidence of the
propagation gap for shear waves at small wave numbers. Later, intensive
studies by molecular dynamics (MD) simulations theory of non-hydrodynamic
shear waves were focused on damping effects and different aspects of
transverse dynamics of liquidsAld84 ; Mac84 ; Bry00 ; Bry11 . The developed
theory Mac84 ; Bry00 ; Bry11  explicitly connected the emergence of shear
waves and their damping with non-hydrodynamic shear-stress relaxation.
Surprisingly, in 2017 Trachenko et al Yan17  obtained from a macroscopic
continuum equation the same expression for the propagation gap of shear waves,
and claimed that ill-defined ”Frenkel jumps” and their characteristic single-
particle ”Frenkel time” are responsible for the propagation gap, that was
immediately criticized inBry18 .

For many-component liquids the list of non-hydrodynamic collective modes is
much larger, because typical high-frequency overdamped excitations for binary
and many-component liquids are optic-like collective modes Bry99 ; Bry00 ;
Bry02 ; Bry11 . Generalized hydrodynamics for the case of transverse dynamics
in binary liquids was developed in Bry00b , and in Bry00  were derived
expressions for the dispersion and damping of transverse optic-like modes in
the long-wavelength limit within the formalism of generalized collective modes
(GCM) Mry95 . It was revealed that the mutual diffusivity and tendency to
demixing (large values of the concentration-concentration structure factor
$S_{xx}(k\to 0)$) define the damping and can cause a similar propagation gap
in $k$-space in case of larger value than the critical damping for optic
modes. And it became quite clear why for ionic liquids the optic modes are
well defined - because of $S_{xx}\stackrel{{\scriptstyle k\to 0}}{{\to}}0$ due
to electroneutrality condition practically vanishing damping of optic modes in
the long-wavelength limit was observed for molten saltsBry04 . The dispersion
curves of longitudinal and transverse excitations for binary liquids were
intensively studied within the GCM methodologyBry00 ; Bry02 ; Bry04 ; Mry09 ,
however, the issue of the width of propagating gap for shear waves in binary
and many-component liquids was not studied analytically. There exists just the
standard propagation gap width Mac84 ; Bry00

$k_{s}=\frac{\sqrt{\rho G_{\infty}}}{2\eta}~{},$ (2)

where $G_{\infty}$ is the macroscopic high-frequency shear modulus. However,
the shear waves as non-hydrodynamic excitations couple with optic modes in
many-component liquids and this coupling should affect the location of
$k_{s}$. In order to obtain the analytical expression for the gap width in
this case one has to solve analytically in the long-wavelength limit a four-
variable dynamic model, which incorporates two two-variable dynamic models
solved in Bry00  separately. The numerical GCM approach for the gap width for
shear waves in binary liquids points out an increase of the gap width with the
mass component ratio, while an increase of concentration of heavier particles
shows the opposite tendencyMry09 . Therefore it will be of great interest to
obtain analytical expressions to explain the simulations and GCM numerical
studies.

The remaining part of the paper is organized as follows. In the next Section
we report our MD simulations of four binary liquids with different mass ratio,
which however have identical static structure, and numerical GCM results for
dispersions oftwo branches of transverse collective excitations. In Section
III we formulate the four-variable dynamic model for transverse dynamics of a
binary liquid, generate a generalized hydrodynamic matrix for it and obtain
four dynamic eigenmodes in the long-wavelength limit as well as an expression
for the gap width $k_{s}$ for shear waves. The last Section contains
conclusion of this study.

##  II Computer simulations of transverse dynamics in a Kob-Andersen liquid
with different mass ratio

The idea of this study is to rationalize an effect of the mass ratio increase
on transverse dynamics in binary liquids, and on the propagation gap for shear
waves in particular. For longitudinal dynamics in binary liquids such an
effect was studied in Bry05 , while for the case of transverse dynamics
Mryglod et al Mry09  reported observation of the increasing gap width for
shear waves vs. mass ratio, although no theory explaining that finding was
developed. Here we will study a binary liquids with identical static structure
but different mass ratio of components $R=m_{B}/m_{A}$, will obtain the
propagation gap dependence on $R$ and develop a theory able to explain the
observed effect.

First we simulated the standard Kob-Andersen 80-20 binary liquid with mass
ratio $R=1$ Kob95 ; Ped18  at temperature $T^{*}=2.0$ and density $n^{*}=1.2$.
The standard set of parameters for Lennard-Jones potential of the KA modes was
used (in relative units to A-A interaction): $\varepsilon_{AA}^{*}=1.0$,
$\varepsilon_{AB}^{*}=1.5$, $\varepsilon_{BB}^{*}=0.5$, $\sigma_{AA}^{*}=1.0$,
$\sigma_{AB}^{*}=0.8$, $\sigma_{BB}^{*}=0.88$. All the output quantities are
in dimensions of SI units. For simulations we have chosen
$\varepsilon_{AA}=119.8~{}K$, $\sigma_{AA}=3.405$ Åand ${\bar{m}}=7.960\cdot
10^{-26}kg$. We have fixed the total mass density of the system and generated
starting configurations with additional 3 mass ratio $R=2$, $R=8$, $R=20$,
then set the total momentum of each starting configuration to zero and
equilibrated each system in microcanonical ensemble till the average
temperature was fluctuating around $T^{*}$. Note, that prepared in this way
the 4 simulated system have identical static structure but different dynamics.
One can see in Fig.1 that there is no difference between three partial pair
distribution functions shown for the case of mass ratio $R=1$ and $R=20$ (for
others, $R=2$ and $R=8$, this is the same).

Figure 1: Identical static structure from pair distribution functions for Kob-
Andersen liquid mixture at density $n^{*}=1.20$ and temperature $T^{*}=2.0$
with different mass ratio: R$=1$ (lines) and R$=20$ (plus symbols).

For each system with corresponding mass ratio we performed production runs of
600 000 time steps and saved time evolution of the following dynamic
variables: spatial Fourier-components of partial transverse mass-currents

${\bf
J}^{T}_{\alpha}(k,t)=\frac{1}{\sqrt{N}}\sum_{j=1}^{N_{\alpha}}m_{\alpha}{\bf
v}^{T}_{j,\alpha}(t)e^{-i{\bf kr}_{j,\alpha}(t)}~{},\alpha=A,B,$ (3)

where ${\bf v}^{T}_{j,\alpha}(t)=[{\bf k}\times{\bf v}_{j,\alpha}]/k$ is the
transverse projection (transverse component) of velocity ${\bf
v}_{j,{\alpha}}$ of $j$-th particle of kind $\alpha$, and their first time
derivatives

$\frac{d{\bf J}^{T}_{\alpha}(k,t)}{dt}\equiv{\dot{\bf
J}}^{T}_{\alpha}(k,t)=\frac{1}{\sqrt{N}}\sum_{j=1}^{N_{\alpha}}[{\bf
F}^{T}_{j,\alpha}(t)-im_{\alpha}({\bf kv}_{j,\alpha}){\bf
v}^{T}_{j,{\alpha}}(t)]e^{-i{\bf kr}_{j,\alpha}(t)}~{},$ (4)

where ${\bf F}^{T}_{j,\alpha}(t)$ is the transverse component of force acting
on the $j$-th particle of kind $\alpha$. The four partial dynamic variables
form a basis set

$A^{(4Tp)}(k,t)=\\{{\bf J}^{T}_{A}(k,t),{\bf J}^{T}_{B}(k,t),{\dot{\bf
J}}^{T}_{A}(k,t),{\dot{\bf J}}^{T}_{B}(k,t)\\}~{},$ (5)

which is used for construction of the $4\times 4$ generalized hydrodynamic
matrix for transverse dynamics of a binary liquid Bry00b . All the matrix
elements (static averages or correlation times) were calculated directly from
MD data. All the $k$-dependent quantities were averaged over all possible
directions of the ${\bf k}$ vectors with the same absolute value.

The transverse total current autocorrelation functions
$F^{T}_{J_{t}J_{t}}(k,t)$ for binary liquids at fixed density and temperature
but with different mass ratio have different relaxation, as it is shown in
Fig.2 at $k=0.302$Å-1. For $R=1$ and $R=2$ the relaxation shape of
$F^{T}_{J_{t}J_{t}}(k,t)$ is typical like in the $k$-region with propagating
shear waves, while for mass ratio $R=8$ and $R=20$ the shape is close to the
hydrodynamic one (1).

Figure 2: Transverse total current autocorrelation function at $k=0.302$Å-1
for Kob-Andersen liquid mixture with different mass ratio R$=m_{B}/m_{A}$ at
$T^{*}=2.0$. The timescale $\tau$ is 3.97299 ps.

This finding stimulated us to calculate the dispersion curves of GCM
propagating eigenmodes for binary liquids at fixed density but with different
mass ratio R using the $4\times 4$ generalized hydrodynamic matrices generated
on the 4-variable set (5). In Fig.3 we show the dispersion of two branches of
transverse propagating modes: low-frequency shear waves and high-frequency
optic-like modes (for which the nearest neighbors move in opposite directions
in transverse planes with respect to the propagation vector ${\bf k}$). The
most striking feature in Fig.3 is the strong effect of the mass ratio on the
low-frequency branch of shear waves: the width of their propagation gap
increases with $R$ and for large mass ratio of components one can expect the
low-frequency branch of shear waves to disappear totally from the frequency
spectrum.

Figure 3: Dispersion of transverse optic and shear wave eigenmodes, obtained
by GCM method, for Kob-Andersen liquid mixture with different mass ratio R at
$T^{*}=2.0$.

In order to check how these findings agree with the expression for the gap
width for shear waves (2) we calculated the shear-stress autocorrelation
functions and estimated shear viscosity for the 4 simulated system with fixed
density but different mass ratio. As expected, the identical static structure
resulted in identical (within error bars) values of the high-frequency shear
modulus $G_{\infty}=3.36\pm 0.02~{}GPa$. The shear viscosity, though, was
obtained from the standard Kubo integration of the shear-stress
autocorrelation functions and showed a small decrease with increasing mass
ratio $R$ (Fig.4). In Bry05  longitudinal collective dynamics was studied for
equimolar binary liquids with fixed mass density and different mass ratio, and
it was observed the increase of mutual diffusivity vs mass ratio $R$. Hence,
it is consistent with the actual observation that the shear viscosity was
decreasing with $R$ as it is shown in Fig.4.

Figure 4: Shear viscosity as a function of mass ratio R, obtained via Kubo
intergrals of shear-stress autocorrelation functios for Kob-Andersen liquid
mixture at $n^{*}=1.2$ and $T^{*}=2.0$.

Now, having the values of $\eta$ and $G_{\infty}$ we can estimate how the
width of propagation gap for shear waves changes according to Eq. (2) and from
MD simulations. In Fig.5 we show the both dependences,and it is seen that the
expression for gap width is not working for binary liquids, especially with
mass ratio larger than 2. Hence, there is a need to develop a theory for shear
waves in binary liquids accounting for the coupling between two high- and low-
frequency branches of collective excitations (i.e. cross-correlations between
total and mass-concentration fluctuations).

Figure 5: Gap width for shear waves vs. mass ratio R for a binary Kob-Andersen
liquid with density $n^{*}=1.2$ at $T^{*}=2.0$ as obtained from Eq.1 (plus
symbols) and observed from MD simulations (cross symbols).

The static cross-correlation between the transverse total and mass-
concentration currents is zero as it was shown in Bry00 . However, the cross-
correlation exists between their first time derivatives, and we can directly
calculate them from MD simulations. In Fig.6 we show the evidence that the
directly estimated from MD cross-correlations
$\langle{\dot{J}}^{T}_{t}(-k){\dot{J}}^{T}_{x}(k)\rangle$ are proportional in
the long-wavelength limit to $k^{2}$, and increase with increasing mass ratio
$R$. Having this information we will try to solve analytically a 4-variable
dynamic model for description of the transverse dynamics in binary liquids in
order to derive the equation for the propagation gap $k_{s}$ of shear waves.

Figure 6: Cross-correlations
$\langle{\dot{J}}^{T}_{t}(-k){\dot{J}}^{T}_{x}(k)\rangle/\langle
J^{T}_{t}(-k)J^{T}_{t}(k)\rangle$ as a function of wave number at different
mass ratio R. Symbols - MD values, fitting solid lines - evidence of the
$k^{2}$-dependence in the long-wavelength region with increasing coupling
constant $a$ vs R.

##  III Four-variable theory of generalized collective modes for transverse
dynamics

We start from the choice of a set of dynamic variables for defining the
generalized hydrodynamic matrix ${\bf T}(k)$ for analytical study of
transverse dynamics of binary liquids. The eigenmodes of ${\bf T}(k)$ allow to
represent the time correlation functions of interest via a weighted sum of
mode contributions in a wide range of wave numbers $k$Mry95 :

$F^{T}_{JJ}(k,t)=\sum_{i=1}^{N_{v}}G_{JJ}^{i}(k)e^{-z_{i}(k)t}~{},$ (6)

where $z_{i}(k)$ are the (complex in general case) eigenvalues of the
$N_{v}\times N_{v}$ matrix ${\bf T}(k)$ for the transverse dynamics, and the
weight coefficients (complex in general case) $G_{JJ}^{i}(k)$ are estimated
via corresponding eigenvectors Mry95 . The simplest choice of dynamic
variables for transverse dynamics would be to start from partial currents of A
and B components (3), as it is quite simply to sample them in MD simulations
and use in numerical GCM study. Note that the GCM eigenvalues will be
identical for matrices ${\bf T}(k)$ generated either on a set of partial
dynamic variables or on a set of any their linear combination. For analytical
study it is better to start from hydrodynamic variables which correspond to
fluctuations of the conserved quantities. Therefore the hydrodynamic variable
of transverse component of total mass current

${\bf J}^{T}_{t}(k,t)=\sum_{\alpha}{\bf J}^{T}_{\alpha}(k,t)~{},\alpha=A,B,$
(7)

is simply the sum of partial currents $J^{T}_{\alpha}(k,t)$, and one can form
from the partial currents another dynamic variable of mass-concentration
current

${\bf J}^{T}_{x}(k,t)=x_{B}{\bf J}^{T}_{A}(k,t)-x_{A}{\bf J}^{T}_{B}(k,t)~{},$
(8)

which is orthogonal to the hydrodynamic variable, i.e.

$\langle{\bf J}_{x}^{T}(-k){\bf J}_{t}^{T}(k)\rangle\equiv 0~{}.$

In (8) $x_{A}=m_{A}N_{A}/(m_{A}N_{A}+m_{B}N_{B})$ and $x_{B}=1-x_{A}$ are
mass-concentrations. Having two orthogonal dynamic variables (7) and (8) we
can extend the inital set of these two ones by other orthogonal dynamic
variables, constructed as it is standard in the GCM approach via first time
derivatives of (7) and (8). Hence, we have obtained a set of four variables

$A^{(4T)}(k,t)=\\{{\bf J}^{T}_{t}(k,t),{\bf J}^{T}_{x}(k,t),{\dot{\bf
J}}^{T}_{t}(k,t),{\dot{\bf J}}^{T}_{x}(k,t)\\}~{},$ (9)

which will be used for generating the generalized hydrodynamic matrix for the
case of transverse dynamics in binary liquids. Note, that inBry00  the
analytical solutions were obtained for two separated (uncoupled) 2-variable
dynamic models: for description of dynamics of total (t) mass-current
fluctuations

$A^{(2Tt)}(k,t)=\\{{\bf J}^{T}_{t}(k,t),{\dot{\bf J}}^{T}_{t}(k,t)\\}$ (10)

and for mass-concentration (x) current fluctuations

$A^{(2Tx)}(k,t)=\\{{\bf J}^{T}_{x}(k,t),{\dot{\bf J}}^{T}_{x}(k,t)\\}~{}.$
(11)

The two-variable dynamic model (10) is equivalent to description of the
transverse dynamics of one-component liquid and results in the long-wavelength
limit in two relaxing (purely real eigenvalues) modes: the standard
hydrodynamic mode (1)

$z_{1}(k)\equiv Re\\{z_{1}(k)\\}=D_{v}k^{2}$ (12)

and a non-hydrodynamic shear-stress relaxation mode Bry11  with the purely
real eigenvalue

$z_{2}(k)\equiv Re\\{z_{2}(k)\\}=d_{0}-D_{v}k^{2}$ (13)

with

$d_{0}=\frac{G_{\infty}}{\eta}=\tau_{M}^{-1}$ (14)

where $\tau_{M}$ is the standard Maxwell relaxation time. The dispersion of
shear waves Chu69  which emerge outside the hydrodynamic regime of pure
liquids is as follows:

$z^{\pm}_{sw}(k)=\sigma_{sw}(k)\pm\sqrt{\langle\omega^{2}_{tt}\rangle(k)-\sigma^{2}_{sw}(k)},$
(15)

with the normalized second frequency momentum of total transverse mass-current
spectral function
$\langle\omega^{2}_{tt}\rangle(k)=\frac{k^{2}G_{\infty}(k)}{\rho}$ and damping
$\sigma_{sw}(k)$, where

$c_{T}^{2}=\frac{G_{\infty}}{\rho}$

is the square of propagation speed for macroscopic ”bare” (non-damped)
transverse excitations. Note, that $c_{T}$ depends on the total mass density
and high-frequency shear modulus, and does not depend explicitly on the ratio
of mass components in the mixture.

Let us apply the set of four dynamic variables (9) to construction of the
$4\times 4$ generalized hydrodynamic matrix ${\bf T}^{(4T)}(k)$ following the
methodology suggested in Mry95 . We result in

$\begin{array}[]{l}{\bf T}^{(4T)}(k)=\left(\begin{array}[]{cccc}0&0&-1&0\\\
0&0&0&-1\\\
c_{T}^{2}k^{2}&ak^{2}&\frac{c^{2}_{T}}{D_{v}}&\frac{ak^{2}D_{12}S_{xx}}{x_{1}x_{2}k_{B}T}\\\
x_{1}x_{2}ak^{2}&\langle\omega_{xx}^{2}\rangle&\frac{x_{1}x_{2}a}{D_{v}}&\frac{\langle\omega_{xx}^{2}\rangle
D_{12}S_{xx}}{x_{1}x_{2}k_{B}T}\\\ \end{array}\right)\ ,\end{array}$ (16)

where the cross-correlations between the first time derivatives of total and
mass-concentration currents was represented as
$\langle\dot{J}_{x}^{T}(-k)\dot{J}_{t}^{T}(k)\rangle=x_{1}x_{2}\bar{m}k_{B}Tak^{2}$
with a coupling constant $a$. Note, that when $a=0$ our $4\times 4$
generalized hydrodynamic matrix decomposes into two $2\times 2$ independent
blocks representing separately dynamics of total transverse mass-current

$\begin{array}[]{l}{\bf T}^{(2Tt)}(k)=\left(\begin{array}[]{cc}0&-1\\\
c_{T}^{2}k^{2}&\frac{c^{2}_{T}}{D_{v}}\\\ \end{array}\right)\ ,\end{array}$
(17)

and mass-concentration transverse current

$\begin{array}[]{l}{\bf T}^{(2Tx)}(k)=\left(\begin{array}[]{cc}0&-1\\\
\langle\omega_{xx}^{2}\rangle&\frac{\langle\omega_{xx}^{2}\rangle
D_{12}S_{xx}}{x_{1}x_{2}k_{B}T}\\\ \end{array}\right)\ ,\end{array}$ (18)

within a precision of the second frequency moment of corresponding transverse
current spectral function. Analytical solutions for eigenmodes of matrices
(17) and (18) were reported in Bry00 . $\langle\omega_{xx}^{2}\rangle$ is the
second frequency moment of mass-concentration transverse current spectral
function.

For the $(4\times 4)$ eigenvalue problem with the generalized hydrodynamic
matrix (16) we are trying to obtain the long-wavelength eigenvalues in a form
of expansion in $k$ restricting the highest order of $O(k^{2})$. Hence we
obtained two purely real eigenvalues

$\begin{array}[]{l}z^{(4T)}_{1}(k)\equiv Re{z^{(4T)}_{1}(k)}=D_{v}k^{2}\\\
z^{(4T)}_{2}(k)\equiv Re{z^{(4T)}_{2}(k)}=d_{0}-D_{v}k^{2}-\Delta
a^{2}k^{2}~{}.\end{array}$ (19)

It is important to note, that the hydrodynamic relaxation mode
$z^{(4T)}_{1}(k)$ is not affected by the non-hydrodynamic processes as it
should be. The other two eigenvalues correspond to a pair of propagating
eigenmodes $z^{(4T)}_{3,4}(k)\equiv z_{opt}^{\pm}(k)$, where

$\begin{array}[]{l}z_{opt}^{\pm}(k)=\sigma_{opt}(k)\pm i\omega_{opt}(k)\\\
\sigma_{opt}(k)=\sigma_{0}+\frac{\Delta}{2}a^{2}k^{2}\\\
\omega_{opt}(k)=\sqrt{\langle\omega_{xx}^{2}\rangle_{k}-\sigma_{opt}^{2}(k)}~{},\end{array}$
(20)

which are the optic-like modes. In (19) and (20) we made a shortcut $\Delta$
for

$\Delta=\frac{x_{1}x_{2}}{D_{v}}\frac{1-\frac{2\sigma_{0}d_{0}}{\langle\omega_{xx}^{2}\rangle_{k=0}}}{(z^{+}_{opt}(k=0)-d_{0})(z^{-}_{opt}(k=0)-d_{0})}~{},$
(21)

and $\sigma_{opt}(k)$ and $\omega_{opt}(k)$ in (20) are the damping and
frequency of the optic modes, which have the following nonzero values in the
long-wavelength limit:

$\begin{array}[]{l}\sigma_{0}=\frac{\langle\omega_{xx}^{2}\rangle_{k=0}D_{12}S_{xx}}{2x_{1}x_{2}k_{B}T}\\\
\omega_{0}=\sqrt{\langle\omega_{xx}^{2}\rangle_{k=0}-\sigma_{0}^{2}}\end{array}$
(22)

It is easily to check the sum rule, that the sum of all four eigenvalues (19)
and (20) equals the sum of the diagonal matrix elements of ${\bf
T}^{(4T)}(k)$. Having the solutions for optic branch (20) we can exclude them
from the characteristic polynomial of the eigenvalue problem for $(4\times 4)$
matrix ${\bf T}^{(4T)}(k)$ (16) and obtain an effective quadratic equation

$z^{2}-(c_{T}^{2}/D_{v}-\Delta
a^{2}k^{2})z-k^{2}(x_{1}x_{2}a^{2}k^{2}-c_{T}^{2}\langle\omega_{xx}^{2}\rangle)/\langle\omega_{xx}^{2}\rangle=0~{},$
(23)

which immediately results in two eigenmodes

$z_{1,2}(k)=\delta(k)\pm\sqrt{\delta^{2}(k)+k^{2}(\frac{x_{1}x_{2}a^{2}k^{2}}{\langle\omega_{xx}^{2}\rangle}-c^{2}_{T})}$
(24)

where $\delta(k)=(\frac{c_{T}^{2}}{D_{v}}-\Delta a^{2}k^{2})/2$. Here the term
with $\Delta$ appears due to coupling of the non-hysrodynamic shear-stress
relaxation with optic modes. It is straightforward to check that when $k\to 0$
one obtains exactly the two relaxing modes (12) and (13). The shear waves
emerge in the binary liquid when the expression under square root becomes
negative and $z_{1,2}(k)$ in (24) become a pair of complex-conjugated numbers.
One can see, that the second term under square root in (24) is a difference of
a positive function of $k^{4}$ and $c^{2}_{T}k^{2}$. The latter is the square
of ”bare” dispersion of non-damped transverse sound-like modes with the
propagation speed $c_{T}$, while the former renormalizes the ”bare” dispersion
down preventing appearence of damped shear waves especially for large values
of $t-x$ cross-correlation coupling $a$. The equation for the gap for shear
waves in binary liquids follows from the expression under square root in (24)
and reads:

$(\frac{x_{1}x_{2}a^{2}}{\langle\omega_{xx}^{2}\rangle}+\frac{\Delta^{2}a^{4}}{4})k_{s}^{4}-(c^{2}_{T}+d_{0}\Delta
a^{2})k_{s}^{2}+\frac{d_{0}}{4}=0~{},$ (25)

which is a biquadratic equation for the propagation gap $k_{s}$. One can
easily make sure that when the $t-x$ cross-correlations are absent, i.e.
$a=0$, one obtains exactly the expression (2) for $k_{s}$. Hence, the effect
of mass-concentration fluctuations on the gap for shear waves is defined by
the $t-x$ coupling constant $a$ and parameter $\Delta$ (21) which accounts for
the difference between the non-hydrodynamic eigenmodes (shear-stress
relaxation and optic modes) in the long-wavelength limit.

##  IV Conclusion

We performed a study of transverse collective excitations in four binary Kob-
Andersen liquids at the same density $n^{*}=1.2$ and temperature $T^{*}=2.0$
but with different mass ratio R. The spectra of transverse collective
excitations in binary liquids contain two branches: high-frequency one of
transverse optic modes and low-frequency one of shear waves. A striking
feature was observed in the behavior of the low-frequency branch: the width of
propagation gap of shear waves increased with the mass ratio R. For very large
mass ratio R there is the tendency of the disappearence of the low-frequency
branch, that is in agreement with previous simulation and GCM results Mry09 .

The standard expression for the gap width for shear waves, derived for pure
liquids is not able to explain the onset of non-hydrodynamic transverse
excitations in binary liquids especially at large mass ratio of components. We
solved analytically a 4-variable dynamic model with accounting for ”t-x”
cross-correlations between total mass and mass-concentration current
correlations. We derived a new equation for the gap width of shear waves in
binary liquids, which reduces to the known gap width expression in case of
absent ”t-x” cross-correlations.

Our findings have another consequence. They show that the dispersion of
collective excitations outside the hydrodynamic regime is strongly dependent
on the mass ratio R, even if the static structure is identical. This raises a
question about the validity of the ”phonon theory of liquid thermodynamics”
Bol12 , which explains the leading contribution to the heat capacity $C_{v}$
as coming from the most short-wavelength longitudinal and transverse non-
damped excitations. According to the ”phonon theory of liquid thermodynamics”
in the case of absent shear waves the heat capacity $C_{v}$ should not be
higher than 2$k_{B}$. However, in our simulations all four binary systems with
identical static structure resulted in the heat capacity $C_{v}=2.70\pm
0.02~{}k_{B}$ that makes evidence of the independence of the heat capacity on
the width of propagation gap of shear waves. Our simulation study of pure and
binary lqiuids on the ”phonon theory of liquid thermodynamics” will be
published elsewhere.

## Acknowledgment

This study was supported by the National Research Foundation of Ukraine under
the grant No. 2020.02/0115.

## References

  * (1) J.-P.Hansen and I.R.McDonald, Theory of Simple Liquids (London: Academic) (1986).
  * (2) J.-P.Boon and S.Yip, Molecular Hydrodynamics (New-York: McGraw-Hill) (1980).
  * (3) C.-H. Chung, S. Yip. Phys. Rev. 182, 323 (1969).
  * (4) B.J. Alder, W.E. Alley. Phys. Today. 37, 56 (1984).
  * (5) R.A. MacPhail, D. Kivelson. J. Chem. Phys. 80, 2102 (1984).
  * (6) T. Bryk, I. Mryglod, J. Phys.:Condens. Matt. 12, 6063 (2000)
  * (7) T. Bryk, Eur. Phys. J. Special Topics 196, 65 (2011); ibid. 227, 2689 (2019).
  * (8) C. Yang, M. T. Dove, V. V. Brazhkin, K. Trachenko, Phys. Rev. Lett. 118, 215502 (2017).
  * (9) T. Bryk, I. Mryglod, G. Ruocco, T. Scopigno, Phys.Rev.Lett. 120, 219601 (2018).
  * (10) T. Bryk, I. Mryglod, Phys.Lett. A. 17, 413 (1999).
  * (11) T. Bryk, I. Mryglod, J. Phys.:Condens. Matt. 14, L445 (2002)
  * (12) T. Bryk, I. Mryglod, Phys. Rev. E. 62, 2188 (2000).
  * (13) I.M. Mryglod, I.P. Omelyan, and M.V. Tokarchuk, Mol. Phys. 84, 235 (1995).
  * (14) T. Bryk, I. Mryglod, J. Phys.:Condens. Matt. 16, L463 (2004)
  * (15) I. M. Mryglod, O. O. Prytula, M. Yu. Druchok. AIP Conference Proceedings 1198, 108 (2009).
  * (16) T. Bryk, I. Mryglod, J. Phys.:Condens. Matt. 17, 413 (2005)
  * (17) W. Kob, H. C. Andersen, Phys. Rev. E 51, 4626 (1995).
  * (18) U.R. Pedersen, T.B. Schrøder, J.C. Dyre, Phys. Rev. Lett. 120, 165501 (2018).
  * (19) D. Bolmatov, V.V. Brazhkin, K. Trachenko, Sci. Rep. 2, 421 (2012).


# Transcending Traditional Boundaries: Leveraging Inter-Annotator Agreement
(IAA) for Enhancing Data Management Operations (DMOps)

Damrin Kim  NamHyeok Kim  Chanjun Park  Harksoo Kim

###### Abstract

This paper presents a novel approach of leveraging Inter-Annotator Agreement
(IAA), traditionally used for assessing labeling consistency, to optimize Data
Management Operations (DMOps). We advocate for the use of IAA in predicting
the labeling quality of individual annotators, leading to cost and time
efficiency in data production. Additionally, our work highlights the potential
of IAA in forecasting document difficulty, thereby boosting the data
construction process’s overall efficiency. This research underscores IAA’s
broader application potential in data-driven research optimization and holds
significant implications for large-scale data projects prioritizing
efficiency, cost reduction, and high-quality data.

Machine Learning, ICML

##  1 Introduction

Within industry settings, data specificity to principal business domains is
crucial, particularly for B2B companies requiring data congruent with customer
needs and business objectivest (Zha et al., 2023). This has led many companies
to create proprietary data and stimulated research into efficient data
construction pipelines, namely Data Management Operations (DMOps) (Choi &
Park, 2023).

In the pursuit of time and cost efficiency in DMOps, preempting labeling
errors to minimize data quality degradation is of paramount importance. This
necessitates effective forecasting of data quality, data production time and
cost, and judicious task allocation amongst workers. Consequently, the ability
to predict data quality in advance to curb unnecessary costs and forestall
worker errors is vital to enhancing data quality (Ge et al., 2017).

In this paper, we present an innovative application of the Inter Annotator
Agreement (IAA) (Artstein, 2017), a metric traditionally used for assessing
labeling consistency, to augment the efficiency of DMOps. We first demonstrate
how IAA can be employed to predict individual workers’ labeling quality. By
measuring and comparing the IAA scores of workers, we identify those whose
scores significantly deviate from the norm, visualized via a heatmap. This
approach facilitates targeted retraining, preemptively mitigating quality
disparities amongst workers and naturally leading to improved data quality,
time, and cost reductions.

Further, we illustrate how IAA can be leveraged to predict document
complexity. By gauging the IAA score of documents during the pilot phase and
comparing these with the scores of pre-existing documents, we can preemptively
forecast document difficulty. This capability enables effective planning for
modeling efforts and task allocation, bolstering the overall operation’s
efficiency.

However, this paper aims to apply IAA not only to measure the level of
agreement in data but also to predict data quality and improve DMOps process
efficiency.

##  2 Predictive Analysis of Labeling Quality by Worker

### Task

In this study, we have chosen to focus on the task of Visual Information
Extraction (VIE) (Yu et al., 2021) underpinned by Optical Character
Recognition (OCR) technology. The complexity of this task necessitates the
labeling of key-value parsing information, enabling the extraction of relevant
data from the recognized text, the understanding of structured information,
and facilitating subsequent analyses or processing.

### Setting

Our experimental design incorporates a variety of official documents. This
includes certificate of resident registration (10 pages), certificate of
archived family relations registration (10 pages), certificate of family
relations (10 pages), and certificate of basic residence registration (10
pages). These documents collectively encompass a wide range of information,
including personal details, family relationships, and specific identification
data such as birth dates, genders, and issuance dates. This information
necessitates tagging by the assigned workers. Our experimental procedure
engaged three workers (A, B, C) for this task.

Figure 1: Kappa coefficients for each worker in the four official documents

### Results

The quality prediction of labeling per worker was examined through the
experimental results delineated in Figure 1. Our primary objective was to
visually render the Inter-Annotator Agreement (IAA) results in a heatmap
format to facilitate an intuitive understanding of discrepancies among various
workers.

As per the heatmap analysis illustrated in Figure 1, there are observable
discrepancies between worker A and workers B & C with respect to certificate
of resident registration and certificate of archived family relations
registration. This observation lends itself to the inference that additional
training and intervention might be requisite for worker A to enhance the
quality of the assigned tasks. Specifically, if the IAA outcomes are
numerically deficient or if significant discrepancies surface among workers,
further work might be necessitated to bolster the reliability and precision of
data quality. This insight could aid in optimizing the efficiency of the DMOps
process by leveraging the IAA.

Conversely, the Kappa coefficient for workers A, B, and C concerning
certificate of family relations and certificate of basic residence
registration exhibits an exceptionally high average of 0.97. This implies that
these datasets are of a very high quality. Hence, we infer that these two
document types may not necessitate any additional work.

In conclusion, IAA serves as a valuable tool for measuring the labeling
quality per worker. If significant inconsistencies among workers are detected
from these measurements, it may be prudent to identify the worker(s) in
question and consider rework. Instead, it should be reserved for situations
where quantitative measurements exhibit substantial discrepancies or fall
below a pre-determined threshold. By adopting this strategic approach to data
quality prediction, we can effectively reduce the time and cost associated
with the work process, leading to a more efficient and streamlined
construction process.

##  3 Analyzing Document Difficulty Prediction

For our second experiment, we compared the IAA between various documents,
including official documents, to determine whether IAA could be used as a
predictor of document difficulty. If document difficulty can be predicted, it
is possible to plan the modeling effort and worker allocation in advance, thus
improving the overall efficiency of the task.

The experiment results are shown in Table 1, with IAA values obtained in the
following order: official documents, diagnoses and medical reports, requests
for quotes, and bills of lading. Generally, documents that are more complex
and contain more parsing information tend to be more difficult to work with,
and the same trend was observed in this experiment. These results suggest that
IAA can not only be used as a measure of consistency but also as a useful
indicator for measuring the difficulty of documents.

Table 1: Comparison and analysis table of IAA results by document

Metrics |  Official
Documents |  Diagnosis and
Medical Report |  Request
for Quote |  Bill
of Lading
---|---|---|---|---
Cohen’s Kappa | 0.976 | 0.967 | 0.927 | 0.772
Fleiss’s Kappa | 0.976 | 0.967 | 0.927 | 0.771
Krippendorff’s Alpha | 0.978 | 0.969 | 0.933 | 0.778

##  4 Conclusions

This paper conducted research to increase the efficiency of DMOps using Inter
Annotator Agreement (IAA). The study demonstrated that IAA can be used for
predicting labeling quality by worker, measuring document difficulty. This
makes it possible to anticipate the difficulty of documents in advance and
plan for modeling efforts and worker allocation accordingly. In other words,
IAA can be used not only as a measure of labeling consistency, but also as a
useful indicator for overall data construction process improvement.

## Acknowledgements

This work was supported by Institute for Information & communications
Technology Planning & Evaluation(IITP) grant funded by the Korea
government(MSIT) (No. 2022-0-00369, (Part 4) Development of AI Technology to
support Expert Decision-making that can Explain the Reasons/Grounds for
Judgment Results based on Expert Knowledge)

## References

  * Artstein (2017) Artstein, R.  Inter-annotator agreement.  _Handbook of linguistic annotation_ , pp. 297–313, 2017.
  * Choi & Park (2023) Choi, E. and Park, C.  Dmops: Data management operation and recipes.  _arXiv preprint arXiv:2301.01228_ , 2023.
  * Ge et al. (2017) Ge, M., O’Brien, T., and Helfert, M.  Predicting data quality success-the bullwhip effect in data quality.  In _Perspectives in Business Informatics Research: 16th International Conference, BIR 2017, Copenhagen, Denmark, August 28–30, 2017, Proceedings 16_ , pp. 157–165. Springer, 2017.
  * Yu et al. (2021) Yu, W., Lu, N., Qi, X., Gong, P., and Xiao, R.  Pick: Processing key information extraction from documents using improved graph learning-convolutional networks.  In _2020 25th International Conference on Pattern Recognition (ICPR)_ , pp. 4363–4370, 2021.  doi: 10.1109/ICPR48806.2021.9412927.
  * Zha et al. (2023) Zha, D., Bhat, Z. P., Lai, K.-H., Yang, F., Jiang, Z., Zhong, S., and Hu, X.  Data-centric artificial intelligence: A survey.  _arXiv preprint arXiv:2303.10158_ , 2023.