\inputencoding{latin1}
\HeaderA{glad}{Analysis of array CGH data}{glad}
\methaliasA{glad.profileCGH}{glad}{glad.profileCGH}
\keyword{models}{glad}
\begin{Description}\relax
This function allows the detection of breakpoints in genomic profiles
obtained by array CGH technology and affects a status (gain, normal
or lost) to each clone.
\end{Description}
\begin{Usage}
\begin{verbatim}
## S3 method for class 'profileCGH':
glad(profileCGH, mediancenter=FALSE,
                smoothfunc="lawsglad", bandwidth=10, round=1.5,
                model="Gaussian", lkern="Exponential", qlambda=0.999,
                base=FALSE, sigma,
                lambdabreak=8, lambdacluster=8, lambdaclusterGen=40,
                type="tricubic", param=c(d=6),
                alpha=0.001, msize=5,
                method="centroid", nmax=8,
                verbose=FALSE, ...)

\end{verbatim}
\end{Usage}
\begin{Arguments}
\begin{ldescription}
\item[\code{profileCGH}] Object of class \code{\LinkA{profileCGH}{profileCGH}}
\item[\code{mediancenter}] If \code{TRUE}, LogRatio are center on their median.
\item[\code{smoothfunc}] Type of algorithm used to smooth \code{LogRatio} by a
piecewise constant function. Choose either \code{lawsglad}, \code{\LinkA{aws}{aws}} or
\code{\LinkA{laws}{laws}}.
\item[\code{bandwidth}] Set the maximal bandwidth \code{hmax} in the
\code{\LinkA{aws}{aws}} or  \code{\LinkA{laws}{laws}} function. For
example, if \code{bandwidth=10} then the \code{hmax} value is set
to 10*\eqn{X_N}{} where \eqn{X_N}{} is the position of the last clone.
\item[\code{round}] The smoothing results are rounded or not depending on
the \code{round} argument. The \code{round} value is passed to the
argument \code{digits} of the \code{\LinkA{round}{round}} function.
\item[\code{model}] Determines the distribution type of the LogRatio. Keep
always the model as "Gaussian" (see \code{\LinkA{laws}{laws}}).
\item[\code{lkern}] Determines the location kernel to be used (see \code{\LinkA{aws}{aws}} or
\code{\LinkA{laws}{laws}}).
\item[\code{qlambda}] Determines the scale parameter for the
stochastic penalty (see \code{\LinkA{aws}{aws}} or
\code{\LinkA{laws}{laws}})
\item[\code{base}] If \code{TRUE}, the position of clone is the physical position onto
the chromosome, otherwise the rank position is used.
\item[\code{sigma}] Value to be passed to either argument \code{sigma2}   
of\code{ \LinkA{aws}{aws}} function or \code{shape} of
\code{\LinkA{laws}{laws}}. If \code{NULL}, sigma is calculated from
the data.
\item[\code{lambdabreak}] Penalty term (\eqn{\lambda'}{}) used during the 
\bold{Optimization of the number of breakpoints} step.
\item[\code{lambdacluster}] Penalty term (\eqn{\lambda*}{}) used during the \bold{MSHR
clustering by chromosome} step.
\item[\code{lambdaclusterGen}] Penalty term (\eqn{\lambda*}{}) used during the \bold{HCSR
clustering throughout the genome} step.
\item[\code{type}] Type of kernel function used in the penalty term during the \bold{Optimization of the
number of breakpoints} step, the \bold{MSHR
clustering by chromosome} step and the \bold{HCSR
clustering throughout the genome} step.
\item[\code{param}] Parameter of kernel used in the penalty term.
\item[\code{alpha}] Risk alpha used for the \bold{Outlier detection} step.
\item[\code{msize}] The outliers MAD are calculated on regions with a
cardinality greater or equal to msize.
\item[\code{method}] The agglomeration method to be used during the \bold{MSHR
clustering by chromosome} and the \bold{HCSR
clustering throughout the genome} clustering steps.
\item[\code{nmax}] Maximum number of clusters (N*max) allowed during
the the \bold{MSHR
clustering by chromosome} and the \bold{HCSR
clustering throughout the genome} clustering steps.
\item[\code{verbose}] If \code{TRUE} some information are printed
\item[\code{...}] 
\end{ldescription}
\end{Arguments}
\begin{Details}\relax
The function \code{glad} implements the methodology which
is described in the article : Analysis of array CGH data: from signal
ratio to gain and loss of DNA regions (Hupé et al., Bioinformatics 2004 20(18):3413-3422).

The principle of the GLAD algorithm:
First, the detection of breakpoints is based on the estimation of a
piecewise constant function with the Adaptive Weights Smoothing (AWS)
procedure (Polzehl and Spokoiny, 2002).
Thus, a procedure based on penalyzed maximum likelihood optimizes the number of
breakpoints allows the undesirable breakpoints to be removed.
Finally, based on the regions previously identified, a two-step
unsupervised classification (\bold{MSHR 
clustering by chromosome} and the \bold{HCSR
clustering throughout the genome}) with model selection criteria
allows a status to be assigned for each region (gain, loss or
normal).

Main parameters to be tuned:
\Tabular{ll}{
\code{qlambda} & if you want the smoothing to fit some very local effect, choose a smaller \code{qlambda}. \\
\code{bandwidth} & choose a bandwidth not to small otherwise you
will have a lot of little discontinuities. \\
\code{lambdabreak} & More the parameter is high more the number
of undesirable breakpoints is high.\\
\code{lambdacluster} & More the parameter is high more the regions
within a chromosome are supposed to belong to the same cluster. \\
\code{lambdaclusterGen} & More the parameter is high more the
regions over the whole genome are supposed to belong to the same cluster.

}
\end{Details}
\begin{Value}
\begin{ldescription}
\item[\code{ }] An object of class "profileCGH" with the following attributes:
\item[\code{profileValues: }] a data.frame with the following added information:

\Itemize{


\item[\bold{Smoothing}] The smoothing values correspond to the median
of each \bold{MSHR (i.e. \code{Region}).}

\item[\bold{Breakpoints}] The last position of a region with identical amount
of DNA is flagged by 1 otherwise it is 0. Note that during the
"Optimization of the number of breakpoints" step, removed
breakpoints are flagged by -1.


\item[\bold{Region}] Each position between two breakpoints are labelled the
same way with an integer value starting from one. The label is
incremented by one when a new breakpoints occurs or when moving to
the next chromosome. The variable \code{region} is what we call MSHR.

\item[\bold{Level}] Each position with equal smoothing value are labelled the
same way with an integer value starting from one. The label is
incremented by one when a new level occurs or when moving to the
next chromosome.

\item[\bold{OutliersAws}] Each AWS outliers are flagged by -1 or 1 
otherwise  it is 0.


\item[\bold{OutliersMad}] Each MAD outliers are flagged by -1 (if it is
in the \eqn{\alpha/2}{} lower tail of the distribution) or 1 (if it is
in the \eqn{\alpha/2}{} upper tail of the distribution)
otherwise  it is 0.

\item[\bold{OutliersTot}] OutliersAws + OutliersMad.

\item[\bold{ZoneChr}] Clusters identified after \bold{MSHR (i.e. \code{Region})
clustering by chromosome}.

\item[\bold{ZoneGen}] Clusters identified after \bold{HCSR clustering throughout the
genome}.

\item[\bold{ZoneGNL}] Status of each clone : Gain is coded by 1, Loss by -1 and
Normal by 0.

}

\item[\code{BkpInfo: }] the data.frame attribute \code{BkpInfo} which gives
the list of breakpoints:
\Itemize{
\item[\bold{PosOrder}] The rank position of each clone on the genome.
\item[\bold{PosBase}] The base position of each clone on the genome.
\item[\bold{Chromosome}] Chromosome name.
}

\item[\code{SigmaC: }] the data.frame attribute \code{SigmaC} gives the estimation of the LogRatio standard-deviation for each chromosome:
\Itemize{
\item[\bold{Chromosome}] Chromosome name.
\item[\bold{Value}] The estimation is based on the Inter Quartile Range.
}

\end{ldescription}
\end{Value}
\begin{Note}\relax
People interested in tools dealing with array CGH analysis can
visit our web-page \url{http://bioinfo.curie.fr}.
\end{Note}
\begin{Author}\relax
Philippe Hupé, \email{glad@curie.fr}.
\end{Author}
\begin{SeeAlso}\relax
\code{\LinkA{profileCGH}{profileCGH}}, \code{\LinkA{as.profileCGH}{as.profileCGH}}, \code{\LinkA{plotProfile}{plotProfile}}.
\end{SeeAlso}
\begin{Examples}
\begin{ExampleCode}

data(snijders)

### Creation of "profileCGH" object
gm13330$Clone <- gm13330$BAC
profileCGH <- as.profileCGH(gm13330)


###########################################################
###
###  glad function as described in Hupé et al. (2004)
###
###########################################################

res <- glad(profileCGH, mediancenter=FALSE,
                smoothfunc="lawsglad", bandwidth=10, round=1.5,
                model="Gaussian", lkern="Exponential", qlambda=0.999,
                base=FALSE,
                lambdabreak=8, lambdacluster=8, lambdaclusterGen=40,
                type="tricubic", param=c(d=6),
                alpha=0.001, msize=5,
                method="centroid", nmax=8,
                verbose=FALSE)

### Genomic profile on the whole genome
plotProfile(res, unit=3, Bkp=TRUE, labels=FALSE, Smoothing="Smoothing",
main="Breakpoints detection: GLAD analysis")

###Genomic profile for chromosome 1
plotProfile(res, unit=3, Bkp=TRUE, labels=TRUE, Chromosome=1,
Smoothing="Smoothing", main="Chromosome 1: GLAD analysis")

### The standard-deviation of LogRatio are:
res$SigmaC

### The list of breakpoints is:
res$BkpInfo

\end{ExampleCode}
\end{Examples}

