source: trunk/doc/stats_weighted.tex @ 480

\documentstyle[12pt]{article}
%\documentstyle[slidesonly]{seminar}

% $Id: error.tex,v 1.11 2004/10/01 09:08:39 peter Exp $

%\input{psfig}

% --- This bit puts ``draft'' over everything! ---
%\special{!userdict begin /bop-hook{gsave 200 30 translate
%65 rotate /Times-Roman findfont 216 scalefont setfont
%0 0 moveto 0.93 setgray (DRAFT) show grestore}def end}
%  --- end of this bit that puts `draft' over everything ---


\flushbottom
\footskip 54pt
\headheight 0pt
\headsep 0pt
\oddsidemargin 0pt
\parindent 0pt
\parskip 2ex
\textheight 230mm
\textwidth 165mm
\topmargin 0pt

\newcommand{\bea} {\begin{eqnarray}}
\newcommand{\eea} {\end{eqnarray}}
\newcommand{\beq} {\begin{equation}}
\newcommand{\eeq} {\end{equation}}
\newcommand{\bibl}[5]
        {#1, {\it #2} {\bf #3} (#4) #5}
\newcommand{\ol}{\overline}

\renewcommand{\baselinestretch} {1.0}
\renewcommand{\textfraction} {0.1}
\renewcommand{\topfraction} {1.0}
\renewcommand{\bottomfraction} {1.0}
\renewcommand{\floatpagefraction} {1.0}

\renewcommand{\d}{{\mathrm{d}}}
\newcommand{\nd}{$^{\mathrm{nd}}$}
\newcommand{\eg}{{\it {e.g.}}}
\newcommand{\ie}{{\it {i.e., }}}
\newcommand{\etal}{{\it {et al.}}}
\newcommand{\eref}[1]{Eq.~(\ref{e:#1})}
\newcommand{\fref}[1]{Fig.~\ref{f:#1}}
\newcommand{\ovr}[2]{\left(\begin{array}{c} #1 \\ #2 \end{array}\right)}

% Use these to include comments and remarks into the text, these will
% obviously appear as footnotes in the final output.
\newcommand{\CR}[1]{\footnote{CR: #1}}
\newcommand{\JH}[1]{\footnote{JH: #1}}
\newcommand{\PE}[1]{\footnote{PE: #1}}
\newcommand{\PJ}[1]{\footnote{PJ: #1}}

\begin{document}

\large
{\bf Weighted Statistics}
\normalsize

\tableofcontents
\clearpage

\section{Introduction}


\section{WeightedAverager}

Building an estimator $T$ of the parameter $\Theta$ there
is a number of criterias that is welcome to be fullfilled.

\begin{itemize}
\item The bias $b=<T> - \Theta$ should be zero. The
estimator is unbiased.

\item The estimator is efficient if the mean squared error
$<(T-\Theta)^2>$ is small. If the estimator is unbiased the
mean squared error is equal to the variance
$<(T-<T>)^2>$ of the estimator.

\item The estimator is consistent if the mean squared error goes to
zero in the limit of infinite number of data points.

\item adding a data point with weight zero should not change the
estimator.

\end{itemize}

We will use these criterias to find the estimator. We will minimize
the variance with the constraint that the estimator must be
unbiased. Also, we check that the estimator is consistent and handles
zero weight in a desired way. 


\subsection{Mean}
We start building an estimator of the mean in the case with varying
variance. 

The likelihood is

\beq
L(m)=\prod (2\pi\sigma_i^2)^{-1/2}e^{-\frac{(x_i-m)^2}{2\sigma_i^2}},
\eeq

which we want to maximize, but can equally well minimize 

\beq
-\ln L(m) = \sum \frac{1}{2}\ln2\pi\sigma_i^2+\frac{(x_i-m)^2}{2\sigma_i^2}.
\eeq

Taking the derivity yields

\beq
\frac{d-\ln L(m)}{dm}=\sum - \frac{x_i-m}{\sigma_i^2}
\eeq

Hence, the Maximum Likelihood method yields the estimator

\beq
m=\frac{\sum x_i/\sigma_i^2}{\sum 1/\sigma_i^2}
\eeq

Let us check the criterias defined above.
First, we want the estimator to be unbiased.

\beq
b=<m>-\mu=\frac{\sum <x_i>/\sigma_i^2}{\sum 1/\sigma_i^2}-\mu=\frac{\sum \mu/\sigma_i^2}{\sum 1/\sigma_i^2}-\mu=0,
\eeq

so the estimator is unbiased.

Second, we examine how efficient the estimator is, \ie how small the
variance is.

\beq
V(m)=\frac{\sum V(x_i)/\sigma_i^2}{\left(\sum 1/\sigma_i^2\right)^2}=\frac{1}{\sum 1/\sigma_i^2},
\eeq

which obviously goes to zero when number of samples (with finite
$\sigma$) goes to infinity. The estimator is consistent.

Trivially we can see that a zero weight data point does not change the
estimator, which was our last condition.

\subsection{Variance}
Let us now examine the case when we do not know the variance, but only
the weight $w_i$ that is proportional to the inverse
of the variance $\sigma_i^2$

\beq
w_i=\frac{\kappa}{\sigma_i^2},
\eeq

and the variance is

\beq
V(m)=\frac{1}{\sum 1/\sigma_i^2}=\frac{\kappa}{\sum w_i}
\eeq

so we need to estimate $\kappa$. The likelihood is now

\beq
L(k)=P(k)\prod \frac{\sqrt{w_i}} {\sqrt{2\pi k}}e^{-\frac{w_i(x_i-m)^2}{2k}},
\eeq

where $P(k)$ is the prior probabilty distribution. If we have no prior knowledge
about $k$, $P(k)$ is constant.

Taking the derivity of the logarithm again yields

\beq
\frac{d-\ln L(k)}{dk}=-\frac{d\ln P(k)}{dk}+\sum \left(\frac{1}{2k}-\frac{w_i(x_i-m)^2}{2k^2}\right)
\eeq

or equivalently

\beq
k=\frac{1}{N}\sum w_i(x_i-m)^2+2k^2\frac{d\ln P(k)}{dk}
\eeq

In principle, any prior probabilty distribution $P(k)$ could be
used. Here we, for simplicity, focus on the one where the last term
becomes a constant, namely

\beq
P(k)=\exp(-\lambda/k).
\eeq

One problem with this choice is that we have to truncate the
distribution in order to normalize it.

The estimator $k$ becomes

\beq
k=\frac{1}{N}\sum w_i(x_i-m)^2+A,
\eeq

where $A$ is constant (depending on $\lambda$ and the truncation point).

Having an estimation of $\kappa$ we can calculate the variance of $m$

\beq
V(m)=\frac{\frac{1}{N}\sum w_i(x_i-m)^2+A}{\sum w_i}=\frac{\frac{1}{N}\sum (x_i-m)^2/\sigma_i^2}{\sum 1/\sigma_i^2}+\frac{A}{\sum 1/\sigma_i^2}
\eeq

Let us now look at estimation of $\kappa$. Is the criterias above
fullfilled? We start looking at the bias

\beq
b=<k>-\kappa=\frac{1}{N}\sum w_i\left<(x_i-m)^2\right>-\kappa
\eeq

Let us look at 

\bea
\left<(x_i-m)^2\right>=
\\\left<(x_i-\mu)^2+(m-\mu)^2-2(x_i-\mu)(m-\mu)\right>
\\V(x_i)+V(m)-2\left<(x_i-\mu)(m-\mu)\right>
\\V(x_i)+V(m)-2\left<(x_i-\mu)(\frac{\sum_j x_j/\sigma_j^2}{\sum_k 1/\sigma_k^2}-\mu)\right>
\\V(x_i)+V(m)-2\left<(\frac{(x_i-\mu)^2/\sigma_j^2}{\sum_k 1/\sigma_k^2})\right>
\\V(x_i)+V(m)-2\frac{1}{\sum_k 1/\sigma_k^2}
\\\sigma_i^2-\frac{1}{\sum_k 1/\sigma_k^2}
\eea

so the bias is

\bea
b=<k>-\kappa=
\\\frac{1}{N}\sum_i (w_i\sigma_i^2-\frac{w_i}{\sum_k 1/\sigma_k^2})-\kappa=
\\\frac{\kappa}{N}\sum_i (1-\frac{w_i}{\sum_k w_k})-\kappa=
\\\frac{\kappa}{N}(N-1)-\kappa
\eea

so the estimator is asymptotically unbiased. If we want the estimation
to be unbiased we could \eg modify $N$ to $N-1$ , exactly as in the
unweighted case, and we get the following estimator of $\kappa$

\beq
k=\frac{1}{N-1}\sum w_i(x_i-m)^2
\eeq

One problem with this estimator is that it is sensitive to weight zero
samples due to the $N$. To solve that we have to express $N$ using
$w$, wich will make the estimator biased. We suggest the substitution

\beq
N\rightarrow\frac{(\sum w_i)^2}{\sum w_i^2}
\eeq

so the estimator finally becomes

\beq
k=\frac{\sum w_i^2}{(\sum w_i)^2-\sum w_i^2}\sum w_i(x_i-m)^2
\eeq

and the variance is

\beq
V(m)=\frac{\sum w_i^2}{(\sum w_i)^2-\sum w_i^2}\frac{\sum w_i(x_i-m)^2}{\sum w_i}+\frac{A}{\sum 1/\sigma_i^2}
\eeq


\section{Score}
\subsection{Pearson}

Pearson correlation is defined as:
\beq
\frac{\sum_i(x_i-\bar{x})(y_i-\bar{y})}{\sqrt{\sum_i (x_i-\bar{x})^2\sum_i (x_i-\bar{x})^2}}.
\eeq
The weighted version should satisfy the following conditions:

\begin{itemize}
\item Having N equal weights the expression reduces to the unweighted case.
\item Adding a pair of data where one the weight is zero does not change the expression.
\item When $x$ and $y$ are identical, the correlation is one. 
\end{itemize}

Therefore we define the weighted correlation to be
\beq
\frac{\sum_iw_i^xw_i^y(x_i-\bar{x})(y_i-\bar{y})}{\sqrt{\sum_iw_i^xw_i^y(x_i-\bar{x})^2\sum_iw_i^xw_i^y(x_i-\bar{x})^2}},
\eeq
where
\beq
\bar{x}=\frac{\sum_i w^x_iw^y_ix_i}{\sum_i w^x_iw^y_i}
\eeq
and
\beq
\bar{y}=\frac{\sum_i w^x_iw^y_iy_i}{\sum_i w^x_iw^y_i}.
\eeq

\subsection{ROC}
If we have a set of values $x^+$ from class + and a set of values
$x^-$ from class -, the ROC curve area is equal to

\beq
\frac{\sum_{\{i,j\}:x^-_i<x^+_j}1}{\sum_{i,j}1}
\eeq

so a natural extension using weights could be

\beq
\frac{\sum_{\{i,j\}:x^-_i<x^+_j}w^-_iw^+_j}{\sum_{i,j}w^-_iw^+_j}
\eeq

\section{Hierarchical clustering}
\label{hc}
A hierarchical clustering consists of two things: finding the two
closest data points, merge these two data points two a new data point
and calculate the new distances from this point to all other points.\\

In the first item, we need a distance matrix, and if we use Euclidean
distanses the natural modification of the expression would be

\beq
d(x,y)=\frac{\sum w_i^xw_j^y(x_i-y_i)^2}{\sum w_i^xw_j^y} \eeq \\

For the second item, inspired by average linkage, we suggest

\beq
d(xy,z)=\frac{\sum w_i^xw_j^z(x_i-z_i)^2+\sum w_i^yw_j^z(y_i-z_i)^2}{\sum w_i^xw_j^z+\sum w_i^yw_j^z}
\eeq

to be the distance between the new merged point $xy$ and $z$, and we
also calculate new weights for this point: $w^{xy}_i=w^x_i+w^y_i$

\section{Regression}
We have the model

\beq 
y_i=\alpha+\beta (x-m_x)+\epsilon_i, 
\eeq 

where $\epsilon_i$ is the noise. The variance of the noise is
inversely proportional to the weight,
$Var(\epsilon_i)=\frac{\sigma^2}{w_i}$. In order to determine the
model parameters, we minimimize the sum of quadratic errors.

\beq
Q_0 = \sum \epsilon_i^2
\eeq

Taking the derivity with respect to $\alpha$ and $\beta$ yields two conditions

\beq
\frac{\partial Q_0}{\partial \alpha} = -2 \sum w_i(y_i - \alpha - \beta (x_i-m_x)=0
\eeq

and

\beq
\frac{\partial Q_0}{\partial \beta} = -2 \sum w_i(x_i-m_x)(y_i-\alpha-\beta(x_i-m_x)=0
\eeq

or equivalently

\beq
\alpha = \frac{\sum w_iy_i}{\sum w_i}=m_y
\eeq

and

\beq
\beta=\frac{\sum w_i(x_i-m_x)(y-m_y)}{\sum w_i(x_i-m_x)^2}
\eeq

Note, by having all weights equal we get back the unweighted
case. Furthermore, we calculate the variance of the estimators of
$\alpha$ and $\beta$.

\beq 
\textrm{Var}(\alpha )=\frac{w_i^2\frac{\sigma^2}{w_i}}{(\sum w_i)^2}=
\frac{\sigma^2}{\sum w_i}
\eeq

and
\beq
\textrm{Var}(\beta )= \frac{w_i^2(x_i-m_x)^2\frac{\sigma^2}{w_i}}
{(\sum w_i(x_i-m_x)^2)^2}=
\frac{\sigma^2}{\sum w_i(x_i-m_x)^2}
\eeq

Finally, we estimate the level of noise, $\sigma^2$. Inspired by the
unweighted estimation

\beq
s^2=\frac{\sum (y_i-\alpha-\beta (x_i-m_x))^2}{n-2}
\eeq

we suggest the following estimator

\beq
s^2=\frac{\sum w_i(y_i-\alpha-\beta (x_i-m_x))^2}{\sum w_i-2\frac{\sum w_i^2}{\sum w_i}}
\eeq

\end{document}