chrBreakpoints {GLAD}R Documentation

Chromosomal breakpoints detection

Description

The function chrBreakpoints.profileCGH find breakpoints corresponding to an abrupt change of the DNA amount along the chromosome.

Usage

chrBreakpoints.profileCGH(profileCGH, smoothfunc="aws", base=TRUE, sigma, bandwidth=10, round=2, verbose=FALSE, ...)

Arguments

profileCGH Object of class "profileCGH".
smoothfunc Type of function used to smooth LogRatio by a piecewise constant function. Choose either aws or laws.
base If TRUE, the position of BAC is the physical position onto the chromosome, otherwise the rank position is used.
sigma Value to be passed to either argument sigma2 of aws function or shape of laws. If NULL, sigma is calculated from the data.
bandwidth Set the maximal bandwidth hmax in the aws or laws function. For example, if bandwidth=10 then the hmax value is set to 10*X_N where X_N is the position of the last BAC.
round The smoothing results of either aws or laws function are rounded or not depending on the round argument. The round value is passed to the argument digits of the round function.
verbose If TRUE some information are printed
... Parameters to be passed to aws or laws function. Typical parameters are qlambda, model, lkern.

Details

The Adaptive Weights Smoothing procedure as described by Polzehl & Spokoiny is used to fit a function piecewise constant. Then, based on the smoothing results of either aws or laws function, a breakpoint is added when two contiguous smoothing values are different: breakpoints are flagged by 1 in the Breakpoints vector and the flag corresponds to the last position of identical amount of DNA. A specific process is implemented for outliers detection: an outlier is a position such that the smoothing value on its right and the smoothing value on its left are equal but the smoothing value in this position is different from the right and left values. If the first position (respectively the last position) is different of the second one (respectively the next to last one) then the position is considered as an outlier. Each outliers are flagged by 1 in the OutliersAws vector.

Value

An object of class "profileCGH" with the following added information in the data.frame attribute profileValues:

Smoothing Smoothing results of either aws or laws function after being rounded or not depending on the round argument.
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.
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.
OutliersAws Each outliers detected are flagged by one otherwise it is 0.
Breakpoints The last position of a region with identical amount of DNA is flagged by 1 otherwise it is 0.

Author(s)

Philippe Hupé, Philippe.Hupe@curie.fr.

References

Polzehl, J. and Spokoiny, S. (2002). Varying coefficient regression modeling by adaptive weights smoothing Manuscript
Polzehl, J. and Spokoiny, S. (2002b). Local likelihood modelling by adaptive weights smoothing, WIAS-Preprint 787
Polzehl, J. and Spokoiny, S. (2000). Adaptive Weights Smoothing with applications to image restoration, J.R.Statist.Soc. B, 62, Part 2, pp. 335-354

See Also

aws, laws.

Examples


data(snijders)
profileCGH <- list(profileValues=gm13330)
class(profileCGH) <- "profileCGH"


# Estimation of the piecewise constant function

res <- chrBreakpoints(profileCGH, smoothfunc="laws",
                      lkern="exponential", model="Gaussian",
                      qlambda=0.999, base=FALSE,bandwidth=10)

plot(LogRatio ~ PosOrder, data=res$profileValues, pch=20)


# Limit between chromosomes

LimitChr <- unique(res$profileValues$LimitChr)+0.5
abline(v=LimitChr, col="grey", lty=2)

lines(res$profileValues$Smoothing ~ res$profileValues$PosOrder, col="green")


# Breakpoints identified

indexBP <- which(res$profileValues$Breakpoints==1)
BP <- res$profileValues$PosOrder[indexBP]+0.5
abline(v=BP, col="red", lty=2)



[Package GLAD version 1.0.1 Index]