\documentclass[xcolor=svgnames]{beamer}
\usepackage[utf8]{inputenc}
\usepackage[english]{babel}
\usepackage{polski}
%\usepackage{amssymb,amsmath}
%\usepackage[latin1]{inputenc}
%\usepackage{amsmath}
%\newcommand\abs[1]{\left|#1\right|}
\usepackage{amsmath}
\newcommand\abs[1]{\left|#1\right|}
\usepackage{hepnicenames}
\usepackage{hepunits}
\usepackage{color}

%% Shortcuts %%                                             
\newcommand{\aver}[1]{\langle #1 \rangle}                   
\newcommand{\dual}[1]{\tilde{#1}}                           
\newcommand{\est}[1]{\widehat{#1}}                          
\newcommand{\ie}{\textit{i.e.}}                             
\newcommand{\nuvec}{\vec{\nu}}                              
\newcommand{\order}[1]{\mathcal{O}\left({#1}\right)}        
\newcommand{\refapp}[1]{app.~(\ref{app:#1})}                
\newcommand{\refeq}[1]{eq.~(\ref{eq:#1})}                   
\newcommand{\reffig}[1]{fig.~(\ref{fig:#1})}                
\newcommand{\refsec}[1]{sec.~(\ref{sec:#1})}                
\newcommand{\reftab}[1]{tab.~(\ref{tab:#1})}                
\newcommand{\rmdx}[1]{\mbox{d} #1 \,} % differential        
\newcommand{\subd}{_{\text{d}}}                             
\newcommand{\subt}{_{\text{t}}}                             
\newcommand{\thvec}{\vec{\vartheta}}                        
\let\oldtheta\theta                                         
\renewcommand{\theta}{\vartheta}                            
\let\eps\varepsilon                                         
\newcommand{\vecest}[1]{\widehat{\vec{#1}}}                 
\newcommand{\wwhat}[1]{\widehat{\widehat{#1}}}              
\newcommand{\xd}{x\subd}                                    
\newcommand{\xt}{x\subt}                                    
\newcommand{\Ekernel}{E(\xd|\xt)}                           
\newcommand{\pdf}{p.d.f.{~}}   

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
\definecolor{mygreen}{cmyk}{0.82,0.11,1,0.25}


\usetheme{Sybila} 

\title[Extracting angular observables with Method of Moments]{Extracting angular observables\\ with Method of Moments}
\author{Marcin Chrz\k{a}szcz$^{1}$ \\in collaboration with \\Frederik Beaujean, Nicola Serra and Danny van Dyk,\\{~}\\ based on \href{http://arxiv.org/abs/1503.04100}{\textit{arXiv:1503.04100}} }
\institute{$^1$~University of Zurich}
\date{\today}

\begin{document}
% --------------------------- SLIDE --------------------------------------------
\frame[plain]{\titlepage}
\author{Marcin Chrz\k{a}szcz{~}}
\institute{(UZH)}
% ------------------------------------------------------------------------------
% --------------------------- SLIDE --------------------------------------------
%\tableofcontents
\begin{frame}
\begin{enumerate}
\item Motivation.
\item Method of Moments.
\item Systematic uncertainties.
\item MC toy studies.
\item Conclusions.
\end{enumerate}

\end{frame}

\placelogotrue
\section{Motivation}
\begin{frame}\frametitle{Motivation}

Likelihood(LL) fits even though widely used suffer from couple of draw backs:
\begin{enumerate}
\item In case of small number events LL fits suffer from convergence problems. This behaviour is well known and was observed several times in toys when we done $\PB \to \PKstar \Pmu \Pmu$.
\item LL can exhibit a bias when underlying physics model is not well known, incomplete or mismodeled.
\item The LL have problems converging when parameters of the \pdf are close to their physical boundaries, so-called ''boundary problem''
\item Accessing uncertainty in LL in some cases requires application of computationally expensive Feldman-Cousins method.
\end{enumerate}



\end{frame}



\begin{frame}\frametitle{Method of Moments}
\begin{columns}
\column{0.05in}{~}
\column{2.2in}
\begin{center} MoM solves the above problems:\end{center}

\column{2.in}
\only<3>{
\begin{center} Drawback:\end{center}
}
\end{columns}



\begin{columns}
\column{0.05in}{~}
\column{2.2in}
\only<1>{
%\begin{center} MoM solves the above problems:\end{center}
\begin{exampleblock}{Advantages of MoM}
  \begin{itemize}
  \item Probability distribution function rapidity converges towards the Gaussian distribution. 
  \item MoM gives an unbias result even with small data sample.
  \item Insensitive to large class of remodelling of physics models.
  \item Is completely insensitive to boundary problems.
  \end{itemize}
\end{exampleblock}
}
\only<2,3>{
%\begin{center} MoM solves the above problems:\end{center}
\begin{exampleblock}{Advantages of MoM}
  \begin{itemize}
  \item Each observable can be determined separately from other.
  \item Uncertainly follows perfectly $1/\sqrt{N}$ scaling.
  \end{itemize}
\end{exampleblock}
}
\column{2.in}
\only<3>{

%\begin{center} Drawback:\end{center}
\begin{alertblock}{Advantages of MoM}
  \begin{itemize}
  \item Estimated uncertainty in MoM is larger then the ones from LL.
  \end{itemize}
\end{alertblock}
}

\end{columns}
\end{frame}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Method of Moments}
\begin{frame}\frametitle{Introduction to MoM}
Let us a define a probability density function \pdf of a decay:
\begin{align}
P(\nuvec, \thvec) \equiv \sum_i S_i(\nuvec) \times f_i(\thvec)
\end{align}
Let's assume further that there exist a dual basis: $\lbrace f_i(\thvec) \rbrace$, $\{\dual{f}_i(\thvec)\}$ that the orthogonality relation is valid:
\begin{equation}                                                                                                                                                                                                                              
    \label{eq:def-ortho-rel}                                                                                                                                                                                                                  
    \int_\Omega \rmdx{\vec{\theta}} \dual{f}_i(\thvec) f_j(\thvec)  = \delta_{ij}                                                                                                                                                             
\end{equation}    
Since we want to use MoM to extract angular observables it's normal to work with Legendre polynomials. In this case we can find self-dual basis: 
\begin{equation}
\forall_i \dual{f}_i = f_i~,
\end{equation}

just by applying the ansatz: $\dual{f}_i=\sum_i a_{ij} f_j$.

\end{frame}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{frame}\frametitle{Determination of angular observables}
Thanks to the orthonormality relation Eq.~\ref{eq:def-ortho-rel} one can calculate the $S_i(\nuvec)$ just by doing the integration:
\begin{align}
S_i(\nuvec)=\int_\Omega d \thvec  P(\nuvec, \thvec) \dual{f}_i(\thvec) 
\end{align}
\pause
We also need to integrate out the $\nuvec$ dependence:
\begin{align}\label{eq:canonical}
\langle S_i \rangle= \int_\Theta d \nuvec  \int_\Omega d \thvec P(\nuvec, \thvec) \dual{f}_i(\thvec) 
\end{align}
\pause
MoM is basically performing integration in~Eq.~\ref{eq:canonical} using MC method:
\begin{align*}
E[S_i] \to \widehat{E[S_i]}=\dfrac{1}{N}\sum_{k=1}^{N} \dual{f}(x_k)
\end{align*}

\end{frame}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{frame}\frametitle{Uncertainty estimation}
MoM provides also a very fast and easy way of estimating the statistical uncertainty:
\begin{align}
\sigma (S_i)= \sqrt{\dfrac{1}{N-1}\sum_{k=1}^N ( \dual{f}_i(x_k) - \widehat{S_i} )^2   } 
\end{align}
and the covariance:
\begin{equation}
\mathrm{Cov}  [S_i, S_j]=\dfrac{1}{N-1} \sum_{k=1}^N [ \widehat{S_i} - \dual{f}_i(x_k) ][ \widehat{S_j} - \dual{f}_j(x_k) ]
\end{equation}
\pause
Thanks to the CLT both equations are satisfied.



\end{frame}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Systematic uncertainties}


%\placelogofalse

\begin{frame}\frametitle{Partial Waves mismodeling}
\begin{columns}
\column{3in}
\only<1>
{

\begin{itemize}
\item Let us consider a decay of $\PB \to P_1 P_2 \Pmuon \APmuon$.
\item In terms of angular \pdf is expressed in terms of partial-wave expansion.
\item For $\PB \to \PK \Ppi \Pmuon \APmuon$ system, S,P,D waves have been studied. 
\end{itemize}

}
\only<2>
{
\begin{itemize}
\item One can write the \pdf separating the hadronic system:
\end{itemize}
\begin{align}
P(\cos \theta_1, & \cos \theta_2, \theta_3)  =  \\  &\quad \nonumber \sum _i S_i(\nuvec, \cos \theta_2) f_i (\cos \theta_1, \theta_3)
\end{align}

}

\column{2in}
\includegraphics[width=0.9\textwidth]{images/fig-topology.pdf}
\end{columns}
\only<1>{
\begin{itemize}
\item The muon system of this kind of decays has a fixed angular dependence in terms of $\theta_1$ and $\theta_3$.
\item The hadron system can have arbitrary large angular momentum.
\end{itemize}
}

\only<2>{

\begin{itemize}
\item $S_i(\nuvec, \cos \theta_2)$ can be further expend in terms of Legendre polynomials $p_l^{\vert m  \vert }(\cos \theta_2)$:
\end{itemize}
\begin{align}
S_i(\nuvec, \cos \theta_2) = \sum_{l=0}^{\inf} S_{k,l}(\nuvec ) p_l^{\vert m \vert}(\cos \theta_2)
\end{align}

\begin{small}


\begin{itemize}
\item Experimentally the $ S_{k,l}$ are easily accessible, but there is a theoretical difficulty as one would need to sum over infinite number of partial waves.
\end{itemize}\end{small}
}
\end{frame}

\placelogotrue

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{frame}\frametitle{Detector effects}
\begin{columns}

\column{2.5in}
\begin{itemize}
\item Since our detectors are not a perfect devices the angular distribution observed by them are not the distributions that the physics model creates.
\end{itemize}
\column{2.5in}
\includegraphics[width=0.9\textwidth]{images/Fig2b.pdf}
\end{columns}
\pause
\begin{itemize}
\item To take into account the acceptance effects one needs to simulate the a large sample of MC events.\\ Try to figure out the efficiency function.
\item Try to figure out the efficiency function.
\item Number of possibilities.
\item Then you can just weight events:
\end{itemize}
\begin{align*}
\widehat{E[S_i]}=\dfrac{1}{\sum_{k=1}^N w_k}\sum_{k=1}^{N} w_k \dual{f}(x_k),~w_k=\dfrac{1}{\epsilon(x_k)}
\end{align*}




\end{frame}

\begin{frame}\frametitle{Unfolding matrix}
\only<1>{
In general one can write the distribution of events after the detector effects:
\begin{align}
P^{\rm{Det}}(x_d) = N \int \int dx_t~ dx_d~ P^{\rm{Phys}}(x_t) E(x_d \vert x_t),
\end{align}
where $N^{-1}=\int \int d x_t~ dx_d~  P^{\rm{Phys}}(x_t) E(x_d \vert x_t)$ and $(x_d \vert x_t)$ denotes the efficiency $\epsilon(x_t)$ and resolution of the detector $ R(x_d\vert x_t)$:
\begin{align}
E(x_d \vert x_t) = \epsilon(x_t) R(x_d\vert x_t)
\end{align}
\pause
One can define the raw moments:
\begin{align}
Q_i^{(m)} = \int \int d x_t~ dx_d~\dual{f}_i(x_d) P^{(m)}(x_t) E(x_d \vert x_t)
\end{align}
The $m$ index corresponds to simulation sample that has $S_0$ and $S_m$ observables set to $\frac{1}{2}$ and rest to zero.
}
\only<2>
{
Once again we can use MC estimator:
\begin{align}
Q_i^{(m)} \to \widehat{Q}_i^{(m)}= \frac{1}{N_t} \sum_i^{N_d} \dual{f}_i(x_d^{i,m})
\end{align}
Linearity of the integral ensures that there has to exists a linear transformation:
\begin{align}
	\label{eq:eqlin}
\vec{Q}=M\vec{S},
\end{align}
where $M$ is so-called unfolding matrix, \pause given by the formula:
\begin{align}
    M_{ij} = \begin{cases}                                                                                                                                   
        2 Q_i^{(0)}                          & j = 0\,,\\                                                                                                    
        2\left(Q_i^{(j)} - Q_i^{(0)}\right)  & j \neq 0\,,\\                                                                                                 
    \end{cases}                                                        
\end{align}
Once we measured the moments $Q$ in data we can invert Eq. 11 and get the $\vec{S}$:
\begin{align}
\widehat{\vec{S}}=M^{-1} \widehat{\vec{Q}},
\end{align}

}
\end{frame}


\section{Toy Studies}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\placelogofalse
\begin{frame}\frametitle{Toy Validation}
\begin{columns}
\column{2.5in}
\begin{itemize}
\item All the statistics properties of MoM have been tested in numbers of TOY MC.
\item As long as you have~$\sim 30$ events your pulls are perfectly gaussian.
\item Uncertainty scales with $\frac{\alpha}{\sqrt{n}}$, $\alpha = \mathcal{O}(1)$.
\item Never observed any boundary problems.
\end{itemize}

\column{2.5in}
\includegraphics[width=0.9\textwidth]{images/pull-Q2_5_6_S5_200.png}\\
\includegraphics[width=0.9\textwidth]{images/Q2_1_2_S5.pdf}

\end{columns}




\end{frame}
\placelogotrue
\begin{frame}\frametitle{Correlation of MoM with Likelihood}

\begin{columns}
\column{2.5in}
\begin{itemize}
\item MoM is highly correlated with LL.
\item Despite the correlation there can be difference of the order of statistical error.
\end{itemize}

\column{2.5in}
\includegraphics[width=0.9\textwidth]{images/S3_scat.pdf}\\


\end{columns}
\begin{columns}
\column{2.5in}
\includegraphics[width=0.9\textwidth]{images/S5_scat.pdf}

\column{2.5in}

\includegraphics[width=0.9\textwidth]{images/S7_scat.pdf}

\end{columns}
\end{frame}

\section{Conclusions}
\placelogotrue
\begin{frame}\frametitle{Conclusions}
\begin{enumerate}
\item MoM posses several big advantages with one drawback which is larger statistical uncertainty.
\item Allows us to go smaller $q^2$ bins (get ready for $1~\GeV^2$ soon!).
\item Alternative method of extracting the detector effects.
\item Can be applied to various rare decays.
\end{enumerate}

\end{frame}

              
\end{document}