normgam.signal: Normal-gamma background correction

In NormalGamma: Normal-gamma convolution model

Description

Performs background correction using the normal-gamma model.

Usage

normgam.signal(x, par, tail.cor = TRUE, cor = 1e-15, gshift = FALSE, 
               mu = par[1], sigma = par[2], k = par[3], theta = par[4])

Arguments

x	vector of observed intensities.
par	vector of parameters; par[1] and par[2] are the mean and standard deviation of the normal distribution and par[3] and par[4] are the shape and scale parameters of the gamma distribution.
tail.cor	logical (see details).
cor	limit of the right tail correction (see details).
gshift	logical; if TRUE and par[3] is smaller than 1, an ad-hoc translation and a thresholding to 0 are applied to background-corrected values so that the mode of corrected value distribution is 0.
mu, sigma	alternative definition of mean and standard deviation of the normal distribution.
k, theta	alternative definition of shape and scale parameters of the gamma distribution.

Details

normgam.signal performs background correction in an additive background noise+signal model with a normal background noise and a gamma-distributed signal. The corrected value from an observed intensity x is the expectation of the signal given the signal and noise distributions. For a set of parameters (mu, sigma, k, theta), it is given by the ratio of the convolution product of dgamma(x, shape=k+1, scale=theta) and dnorm(x, mean=mu, sd=sigma) and the convolution product of dgamma(x, shape=k, scale=theta) and dnorm(x, mean=mu, sd=sigma). For more details see Plancade S., Rozenholc Y. and Lund E., BMC Bioinfo 2012.

If tail.cor = TRUE, a linear approximation of right tail is applied to values with density estimate smaller than cor in the computation of normal-gamma convoluted densities (see dnormgam).

Only one definition of the parameters is required, either par or (mu, sigma, k, theta). If both are specified and do not match, an error message is returned.

Value

Vector of background noise-corrected intensities.

Author(s)

Plancade S. and Rozenholc Y.

References

Plancade S., Rozenholc Y. and Lund E. "Generalization of the normal-exponential model : exploration of a more accurate parametrisation for the signal distribution on Illumina BeadArrays", BMC Bioinfo 2012, 13(329).

See Also

dnormgam computes the density of the normal-gamma distribution and normgam.fit computes the Maximum Likelihood Estimator. Intensities provides an example of Illumina microarray data.

Examples

#Example 1: simulated data

n = 50000
par = c(60,5,0.15,400)
S = rgamma(n, shape=par[3], scale=par[4])
B = rnorm(n, mean=par[1], sd=par[2]) 
X = S + B

par(mfrow=c(2,1))

Shat1 = normgam.signal(X, par)
H1 = histogram(Shat1, type='irregular', verbose=FALSE, plot=FALSE)
plot(H1, xlim=c(0,50))
I = seq(from=0, to=50, length=1000)
lines(I, dgamma(I, shape=0.15, scale=400), col='red')

Shat2 = normgam.signal(X, par, gshift = TRUE)
H2 = hist(Shat2, 10000, plot=FALSE)
plot(H2, xlim=c(0,50), freq=FALSE)
lines(I, dgamma(I, shape=0.15, scale=400), col='red')


#Example 2: Illumina data

## Not run: 

data(RegNegIntensities_Example)

X = Intensities$Regular
N = Intensities$Negative

# parameter estimation
parmle = normgam.fit(X, N)$par 

Shat = normgam.signal(X,parmle)
H = histogram(Shat, type='irregular', verbose=FALSE, plot=FALSE)
plot(H, xlim=c(0,30))

## End(Not run)

NormalGamma documentation built on May 2, 2019, 10:25 a.m.