%\VignetteIndexEntry{STK++: Statistics and Algebra Reference Guide} \documentclass[a4paper,10pt]{article} %------------------------- % preamble for nice lsitings and notes \include{rtkpp-preamble} \geometry{top=2cm, bottom=2cm, left=2cm}%, right=1cm} %% need no \usepackage{Sweave.sty} <>= library(rtkore) rtkore.version <- packageDescription("rtkore")$Version rtkore.date <- packageDescription("rtkore")$Date @ % Title Page \title{Statistics and Algebra Reference Guide} \author{Serge Iovleff} \date{\today} % start documentation \begin{document} \SweaveOpts{concordance=TRUE} \maketitle \begin{abstract} This reference guide gives a general review of the capabilities offered by the \stkpp{} library. The library is divided in various \emph{projects}. The "Arrays" project is described in detail in the vignette "STK++ Arrays, User Guide" (\cite{rtkore-Arrays}) and the quick reference guide. This vignette focus on the "STatistiK" project (which provides statistical tools) and the "Algebra" project (which provide, mainly, matrix decomposition/inversion tools). \end{abstract} \tableofcontents \section{Statistical functors, methods and functions (STatistiK project)} \label{sec:STatistiK} This section describe the main features provided by the STatistiK project. Mainly \begin{enumerate} \item the probability classes (\ref{subsec:prob}), \item the descriptive statistical methods, \item the utilities related methods. \end{enumerate} The computation of factors using as input a vector or a matrix is detailed in section (\ref{sec:factors}). \subsection{Probabilities } \label{subsec:prob} All the probabilities handled by R are available in \rtkore{}. In the stand-alone \stkpp{} library, only a subset of theses probabilities are implemented. Probability distribution classes are defined in the \code{namespace Law} and can be used as in this example \begin{minipage}[t]{0.66\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoProbability.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.33\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoProbability.out} \end{minipage} All probability distribution classes have a similar prototype like the one given below \begin{lstlisting}[style=customcpp,caption=Prototype of probability distribution class (example taken from Cauchy class)] class Cauchy: public IUnivLaw { public: Cauchy( Real const& mu=0, Real const& scale=1); virtual ~Cauchy() {} Real const& mu() const; Real const& scale() const; void setMu( Real const& mu); void setScale( Real const& scale); virtual Real rand() const; // generate a Cauchy random variate virtual Real pdf( Real const& x) const; // probability distribution function (pdf) virtual Real lpdf( Real const& x) const; // log-pdf virtual Real cdf( Real const& t) const; // cumulative distribution function (cdf) virtual Real cdfc( Real const& t) const; // complementary cdf virtual Real lcdf( Real const& t) const; // log-cdf virtual Real lcdfc( Real const& t) const;// log-complementary cdf virtual Real icdf( Real const& p) const; // inverse cdf (quantiles) static Real rand( Real const& mu, Real const& scale); static Real pdf( Real const& x, Real const& mu, Real const& scale); static Real lpdf( Real const& x, Real const& mu, Real const& scale); static Real cdf( Real const& t, Real const& mu, Real const& scale); static Real icdf( Real const& p, Real const& mu, Real const& scale); protected: Real mu_; Real scale_; } \end{lstlisting} If $f$ denote the density of some probability distribution function (pdf) on $\R$, the methods have the following meaning \begin{enumerate} \item \code{pdf(x)} return the pdf value $f(x)$, \item \code{lpdf(x)} return the $\log$-pdf value $\log(f(x))$, \item \code{rand()} return a random variate with pdf $f$, \item \code{cdf(t)} return the cumulative distribution function (cdf) value $F(t)=\int_{-\infty}^t f(x) dx$ \item \code{lcdf(t)} return the $\log$-cdf value $\log F(t)$ \item \code{cdfc(t)} return the complementary cdf value $G(t)=\int_t^{+\infty} f(x) dx$ \item \code{cdfc(t)} return the $\log$-complementary cdf value $\log G(t)$ \item \code{icdf(p)} return the inverse cumulative distribution function value $F^{-1}(p)$. \end{enumerate} The table \ref{tab:prob} gives the list of available probability distribution. \begin{table}[htb] \begin{tabular}{|l|l|l|l|} \hline Name & Constructor & R functions &Notes\\ \hline Bernoulli & \scode{Law::Bernouilli(p) } & - & \\ \hline Binomial & \scode{Law::Binomial(n,p)} & \scode{*binom} & \\ \hline Beta & \scode{Law::Beta(alpha,beta)} & \scode{*beta} & \\ \hline Categorical & \scode{Law::Categorical(p)} & - & p can be any \stkpp{} vector \\ \hline Cauchy & \scode{Law::Cauchy(m,s)} & \scode{*cauchy}& \\ \hline ChiSquared & \scode{Law::ChiSquared(n)} & \scode{*chisq} & \\ \hline Exponential & \scode{Law::Exponential(lambda)} & \scode{*exp} & Parameterization $\lambda e^{-\lambda x}$ \\ \hline FisherSnedecor & \scode{Law::FisherSnedecor(df1,df2)} & \scode{*f} & \\ \hline Gamma & \scode{Law::Gamma(a,b)} & \scode{*gamma} & Parameterization $\frac{x^{a-1}}{\beta^a\Gamma(a)} e^{-x/\beta}$ \\ \hline Geometric & \scode{Law::Geometric(p)} & \scode{*geom} & \\ \hline HyperGeometric & \scode{Law::HyperGeometric(m,n,k)} & \scode{*hyper} & \\ \hline Logistic & \scode{Law::Logistic(mu,scale)} & \scode{*logis} & \\ \hline LogNormal & \scode{Law::LogNormal(mulog,sdlog)} & \scode{*lnorm} & \\ \hline NegativeBinomial & \scode{Law::NegativeBinomial(size,prob,mu)}& \scode{*nbinom}&\\ \hline Normal & \scode{Law::Normal(mu,sigma)} & \scode{*norm} &\\ \hline Poisson & \scode{Law::Poisson(lambda)} & \scode{*poiss} & \\ \hline Student & \scode{Law::Student(df)} & \scode{*t} &\\ \hline Uniform & \scode{Law::Uniform(a,b)} & \scode{*unif} & \\ \hline UniformDiscrete & \scode{Law::UniformDiscrete(a,b)} & - & \\ \hline Weibull & \scode{Law::Weibull(a)} & \scode{*weibull} & \\ \hline \end{tabular} \caption{List of the available probability distribution in \rtkore{}} \label{tab:prob} \end{table} All Distribution laws methods can be applied to a vector/array/expression using the corresponding methods. An example is given below. \begin{lstlisting}[style=customcpp,caption=Compute the log-complementary cdf in a logistic regression] // logis is logistic distribution STK::Law::Logistic logis; // y is a (R) matrix of size (n,p) and beta is a (R) vector of size p STK::RMatrix y; STK::RVector beta; // compute log-complementary cdf (y_*beta).lcdfc(logis); \end{lstlisting} \subsection{Statistical Methods and global functions} \stkpp{} provides a lot of methods, functions and functors in order to compute usual statistics. \subsubsection{Methods} Let $m$ be any kind of array (square, vector, point, etc...). it is possible to compute the $\min$, $\max$, mean, variance of the elements. These computations can be safe (i.e. discarding N.A. and infinite values) or unsafe and weighted \begin{table}[H] \begin{tabular}{|l|l|l|l|} \hline Method & weigthed version & safe versions & Notes\\ \hline \scode{m.norm()} & \scode{m.wnorm(w)} & \scode{m.normSafe(); m.wnormSafe(w)} & $\sqrt{\sum |m_{ij}|^2}$ \\ \hline \scode{m.norm2()} & \scode{m.wnorm2(w)} & \scode{m.norm2Safe(); m.wnorm2Safe(w)} & $\sum |m_{ij}|^2$ \\ \hline \scode{m.normInf()} & \scode{m.wnormInf(w)} & \scode{m.normInfSafe(); m.wnormInfSafe(w)} & $\sup |m_{ij}|$\\ \hline \scode{m.sum()} & \scode{m.wsum(w)} & \scode{m.sumSafe(); m.wsumSafe(w)} & $\sum m_{ij}$ \\ \hline \scode{m.mean()} & \scode{m.wmean(w)} & \scode{m.meanSafe(); m.wmeanSafe(w)} & $\frac{1}{n}\sum m_{ij}$\\ \hline \scode{m.variance()} & \scode{m.wvariance(w)} & \scode{m.varianceSafe(); m.wvarianceSafe(w)} & $\frac{1}{n}\sum (m_{ij}-\bar{m})^2$\\ \hline \scode{m.variance(mu)} & \scode{m.wvariance(mu,w)}& \scode{m.varianceSafe(mu); m.wvarianceSafe(mu,w)}& $\frac{1}{n}\sum (m_{ij}-\mu)^2$\\ \hline \end{tabular} \caption{List of the available statistical methods for the arrays. $n$ represents the number of elements of $m$. For safe versions, $n$ represents the number of available observations in $m$.} \end{table} \subsubsection{Statistical functions} For two dimensional arrays, there exists global functions allowing to compute the usual statistics by column or by row. By default all global functions are computing the statistics columns by columns. For example, if $m$ is an array, \code{sum(m)} return a row-vector of range \code{m.cols()} containing the sum of each columns. The alias \code{sumByCol(m)} can also be used. The sum of each rows can be obtained using the function \code{sumByow(m)} which return an array of range \code{m.rows()}. These computations can be safe (i.e. discarding N.A. and infinite values) or unsafe and/or weighted. \begin{table}[H] \begin{tabular}{|l|l|l|l|} \hline Function & weigthed version & safe versions (\scode{/* */} for optional arg) \\ \hline \scode{min(m)} & \scode{min(m, w)} &\scode{minSafe(m/*,w*/)} \\ \scode{minByCol(m)} & \scode{minByCol(m, w)} &\scode{minSafeByCol(m/*,w*/)} \\ \scode{minByRow(m)} & \scode{minByRow(m, w)} &\scode{minSafeByRow(m/*,w*/)} \\ \hline \scode{max(m)} & \scode{max(m, w)} &\scode{maxSafe(m/*,w*/)} \\ \scode{maxByCol(m)} & \scode{maxByCol(m, w)} &\scode{maxSafeByCol(m/*,w*/)} \\ \scode{maxByRow(m)} & \scode{maxByRow(m, w)} &\scode{maxSafeByRow(m/*,w*/)} \\ \hline \scode{sum(m)} & \scode{sum(m, w)} &\scode{sumSafe(m/*,w*/)} \\ \scode{sumByCol(m)} & \scode{sumByCol(m, w)} &\scode{sumSafeByCol(m/*,w*/)} \\ \scode{sumByRow(m)} & \scode{sumByRow(m, w)} &\scode{sumSafeByRow(m/*,w*/)} \\ \hline \scode{mean(m)} & \scode{mean(m, w)} &\scode{meanSafe(m/*,w*/)} \\ \scode{meanByCol(m)} & \scode{meanByCol(m, w)} &\scode{meanSafeByCol(m/*,w*/)} \\ \scode{meanByRow(m)} & \scode{meanByRow(m, w)} &\scode{meanSafeByRow(m/*,w*/)} \\ \hline \scode{variance(m, unbiased)} & \scode{variance(m, w, unbiased)} &\scode{varianceSafe(m/*,w*/, unbiased)} \\ \scode{varianceByCol(m, unbiased)} & \scode{varianceByCol(m, w, unbiased)} &\scode{varianceSafeByCol(m/*,w*/, unbiased)} \\ \scode{varianceByRow(m, unbiased)} & \scode{varianceByRow(m, w, unbiased)} &\scode{varianceSafeByRow(m/*,w*/, unbiased)} \\ \hline \scode{varianceWithFixedMean(m, mu, unbiased)} & \scode{variance*(m, w, mu, unbiased)} &\scode{variance*Safe(m/*,w*/, mu, unbiased)} \\ \scode{varianceWithFixedMeanByCol(m, mu, unbiased)} & \scode{variance*ByCol(m, w, mu, unbiased)} &\scode{variance*SafeByCol(m/*,w*/, mu, unbiased)} \\ \scode{varianceWithFixedMeanByRow(m, mu, unbiased)} & \scode{variance*ByRow(m, w, mu, unbiased)} & \scode{variance*SafeByRow(m/*,w*/, mu, unbiased)} \\ \hline \end{tabular} \caption{List of the available global statistical functions for arrays. $m$ is the array, $w$ the vector of weights. \code{unbiased} is a Boolean to set to \code{true} if unbiased variance (divided by $n-1$) is desired. } \end{table} The covariance can be computed in two ways : using two vectors of same size, or using all the columns (rows) of a two-dimensional array. In the first case the functions return the value of the covariance, in the second case, they return a \code{CSquareArray}. \begin{table}[H] \begin{tabular}{|l|l|l|l|} \hline Function & weigthed version & safe versions \\ \hline \scode{covariance(v1, v2, unbiased)} & \scode{covariance(v1, v2, w, unbiased)} & \vcell{ \scode{covarianceSafe(v1, v2, unbiased)} \\ \scode{covarianceSafe(v1, v2, w, unbiased)} } \\ \hline \scode{covarianceWithFixedMean(v1, v2, mean, unbiased)} & \scode{covariance*(v1, v2, w, mean, unbiased)} &\scode{covariance*Safe(v1, v2, unbiased)} \\ \hline \hline \scode{covariance(m, unbiased)} & \scode{covariance(m, w, unbiased)} & \\ \scode{covarianceByRow(m, unbiased)} & \scode{covarianceByRow(m, w, unbiased)} & \\ \hline \scode{covarianceWithFixedMean(m, mean, unbiased)} & \scode{covariance*(m, w, mean, unbiased)} & \\ \scode{covarianceWithFixedMeanByRow(m, mean, unbiased)} & \scode{covariance*ByRow(m, w, mean, unbiased)} & \\ \hline \end{tabular} \caption{List of the available covariance functions for vectors and arrays. $v_1$ and $v_2$ are vectors, $m$ is an array, $w$ a vector of weights. \code{unbiased} is a Boolean to set to \code{true} if unbiased covariance (divided by $n-1$) is desired. The first set of covariance functions return a scalar, the second set of covariance functions return a \code{CSquareArray}} \end{table} The following example illustrate the use of the covariance function: \begin{minipage}[t]{0.66\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoStatCovariance.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.33\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoStatCovariance.out} \end{minipage} \subsection{Miscellaneous statistical functions} Given an array $m$, it is possible to center it, to standardize it and to perform the reverse operations. They are listed in table (\ref{tab:MiscStat}) given below. \begin{table}[htb] \begin{tabular}{|l|l|} \hline Function & weighted version \\ \hline \scode{center(m, mean)} & \scode{center(m, w, mean)} \\ \scode{centerByCol(m, mean)} & \scode{centerByCol(m, w, mean)} \\ \scode{centerByRow(m, mean)} & \scode{centerByRow(m, w, mean)} \\ \hline \scode{standardize(m, std, unbiased)} & \scode{standardize(m, w, std, unbiased)} \\ \scode{standardizeByCol(m, std, unbiased)} & \scode{standardizeByCol(m, w, std, unbiased)} \\ \scode{standardizeByRow(m, std, unbiased)} & \scode{standardizeByRow(m, w, std, unbiased)} \\ \hline \scode{standardize(m, mean, std, unbiased)} & \scode{standardize(m, w, mean, std, unbiased)} \\ \scode{standardizeByCol(m, mean, std, unbiased)} & \scode{standardizeByCol(m, w, mean, std, unbiased)} \\ \scode{standardizeByRow(m, mean, std, unbiased)} & \scode{standardizeByRow(m, w, mean, std, unbiased)} \\ \hline \scode{uncenter(m, mean)} & \\ \scode{uncenterByCol(m, mean)} & \\ \scode{uncenterByRow(m, mean)} & \\ \hline \scode{unstandardize(m, std)} & \\ \scode{unstandardizeByCol(m, std)} & \\ \scode{unstandardizeByRow(m, std)} & \\ \hline \scode{unstandardize(m, mean, std)} & \\ \scode{unstandardizeByCol(m, mean, std)} & \\ \scode{unstandardizeByRow(m, mean, std)} & \\ \hline \end{tabular} \caption{List of the available utilities functions for centering and/or standardize an array. $m$ is an array, $w$ a vector of weights. If used on the columns, \code{mean} and \code{std} have to be points (row-vectors), if used by rows, \code{mean} and \code{std} have to be vectors.\code{unbiased} is a Boolean to set to \code{true} if unbiased covariance (divided by $n-1$) is desired.} \label{tab:MiscStat} \end{table} These methods are illustrated in the following example \begin{minipage}[t]{0.66\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoStatTransform.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.33\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoStatTransform.out} \end{minipage} \begin{note} By default, all functions applied on an array are applied column by column. \end{note} \section{Computing factors} \label{sec:factors} Given a finite collection of object in a vector or an array/expression, it is possible to encode it as factor using the classes \code{Stat::Factor} (for vectors) and \code{Stat::MultiFactor} (for arrays). These classes are runners and you have to use the \code{run} method in order to trigger the computation of the factors. An example is given below \begin{minipage}[t]{0.72\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoStatFactors.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.27\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoStatFactors.out} \end{minipage} \section{Linear Algebra classes, methods and functions (Algebra project)} \label{sec:Algebra} \stkpp{} basic linear operation as product, dot product, sum, multiplication by a scalar,... are encoded in template expressions and optimized. Since \stkpp{} version 0.9 and later, \ttcode{\lapack{}} library can be used as back-ends for dense matrix matrix decomposition (QR, Svd, eigenvalues) and least square regression. In order to use \lapack{}, you must activate its usage by defining the following macros \ttcode{-DSTKUSE\lapack{}} at compilation time and by linking your code with your installed \lapack{} library using \ttcode{-l\lapack{}} (at least in *nix operating systems). \begin{table}[H] \begin{tabular}{|l|l|l|} \hline Class & constructor & Note \\ \hline lapack::Qr & Qr( data, ref = false) & if data is an \ttcode{ArrayXX} and \ttcode{ref} is true, \\ (or Qr) & Qr( data) & data will be overwritten by $Q$ \\ \hline lapack::Svd & Svd( data, ref = false, withU=true, withV=true) & if \ttcode{ref} is true, \\ (or Svd) & Svd(data) & data will be overwritten by $Q$ \\ \hline lapack::SymEigen & SymEigen( data, ref = false) & if data is a \ttcode{SquareArray} and \ttcode{ref} is true, \\ (or SymEigen) & SymEigen( data) & data will be overwritten by $Q$ \\ \hline lapack::MultiLeastSquare & MultiLeastSquare( b, a, isBref = false, isAref=false) & \\ (or MultiLeastSquare) & MultiLeastSquare( b, a) & \\ \hline \end{tabular} \end{table} \subsection{Matrix decomposition} All methods for matrix decomposition are enclosed in classes. Computations are launched using the run method. \subsubsection{QR decomposition} $QR$ decomposition (also called a $QR$ factorization) of a matrix is a decomposition of a matrix $A$ into a product $A = QR$ of an orthogonal matrix $Q$ and an upper triangular matrix $R$. $QR$ decomposition can be achieved using either the class \code{STK::Qr} or if \lapack{} is available \code{STK::lapack::Qr}. In later case, code have to be compiled using \ttcode{-DSTKUSE\lapack{}} flag and linked using \ttcode{-llapack} (at least for GNU-like compilers). \begin{minipage}[t]{0.60\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoQr.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.39\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoQr.out} \end{minipage} \begin{note} By default the matrix $Q$ is represented as a product of elementary reflectors $ Q = H_1 H_2 . . . H_k$, where $k = \min(m,n)$ each $H_i$ has the form $H_i = I - \tau vv'$. It is possible to get the $Q$ matrix (of size $(m,m)$) by using the \ttcode{compQ()} method. $Q$ will be overwritten. \end{note} It is possible to update a QR decomposition when a column is added or removed to the original matrix. In the following example, we remove the column number 2 of a matrix and then insert a column with value 1 after the column number 1. \begin{minipage}[t]{0.60\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoQrUpdate.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.39\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoQrUpdate.out} \end{minipage} \subsubsection{SVD decomposition} The singular-value decomposition of an $(m,n)$ real (or complex) matrix $M$ is a factorization of the form $ \mathbf{U\Sigma V^*}$, where $U$ is an $(m,n)$ real (or complex) unitary matrix, $\mathbf{\Sigma}$ is an $(m,n)$ rectangular diagonal matrix with non-negative real numbers on the diagonal, and $\mathbf{V}$ is an $(n,n)$ real (or complex) unitary matrix. The diagonal entries $\sigma_i$ of $\mathbf{\Sigma}$ are known as the singular values of $M$. \begin{minipage}[t]{0.66\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoSvd.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.33\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoSvd.out} \end{minipage} \begin{note} Singular values are stored in a vector in \lapack{} method and in a diagonal matrix in \stkpp{} method. It is also possible to compute only $U$ and/or $V$ matrix. \end{note} \subsubsection{Eigenvalues decomposition} Let $A$ be a square $(n,n)$ matrix with $n$ linearly independent eigenvectors, $ q_i \,\, (i = 1, \dots, N).$ Then $A$ can be factorized as \[ \mathbf{A}=\mathbf{Q}\mathbf{\Lambda}\mathbf{Q}^{-1} \] where $Q$ is the square $(n,n)$ matrix whose $i$-th column is the eigenvector $ q_i$ of $A$ and $\Lambda$ is the diagonal matrix whose diagonal elements are the corresponding eigenvalues, i.e., $ \Lambda_{ii}=\lambda_i $. If $A$ is a symmetric matrix then $Q$ is an orthogonal matrix. STK provide native and \lapack{} interface classes allowing to compute the eigenvalue decomposition of a \emph{symmetric} square matrix. \begin{minipage}[t]{0.66\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoSymEigen.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.33\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoSymEigen.out} \end{minipage} \begin{note} \stkpp{} eigenvalues computation need a full symmetric matrix as input while \lapack{} version use only upper part of the input data. It is also possible to use the lower part of the matrix in the \lapack{} version. \end{note} \subsection{Solving least square problems} In linear regression, the observations are assumed to be the result of random deviations from an underlying relationship between the dependent variables $y$ and independent variable $x$. Given a data set \[ \{y_{i1},\ldots,y_{id},\, x_{i1}, \ldots, x_{ip}\}_{i=1}^n \] of $n$ statistical units, a linear regression model assumes that the relationship between the dependent variable $y$ and the regressors $x$ is linear. This relationship is modeled through a disturbance term that adds "noise" to the linear relationship between the dependent variable and regressors. Thus the model takes the form \[ \mathbf{y}_i = \mathbf{x}^\mathsf{T}_i \boldsymbol\beta + \varepsilon_i, \qquad i = 1, \ldots, n. \] Often these $n$ equations are stacked together and written in matrix notation as \[ \mathbf{Y} = X\boldsymbol\beta + \boldsymbol\varepsilon, \, \] \stkpp{} provide native and \lapack{} interface classes allowing to solve the least square regression problem. \begin{minipage}[t]{0.66\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoLeastSquare.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.33\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoLeastSquare.out} \end{minipage} \begin{note} It is also possible to solve weighted least-square regression problems. \end{note} \subsection{Inverting matrices} Matrices can be inverted using either the templated functor STK::InvertMatrix or the templated function STK::invert. The first example below deals with general square and/or symmetric matrices. \begin{minipage}[t]{0.55\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoInvert.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.43\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoInvert.out} \end{minipage} The second example deals with lower and upper triangular arrays. \begin{minipage}[t]{0.66\textwidth} \lstinputlisting[style=customcpp,caption=Example]{programs/tutoInvertTriangular.cpp} \end{minipage} \hspace{0.2cm} \begin{minipage}[t]{0.33\textwidth} \addtocounter{lstlisting}{-1} \lstinputlisting[style=customcpp,caption=Output]{programs/tutoInvertTriangular.out} \end{minipage} \begin{note} If the matrix is not inversible, the result provided will be a generalized inverse. \end{note} %----------------------------------------- \bibliographystyle{plain} \bibliography{rtkore} \end{document}