The code is well commented for further information.
First, genesimulator
is called to obtain a vector of
mean gene intensities (for a number of genes and a number of replicates
for each gene.
Then link{simdurbin2}
simulates a series of gene intensities
using the (log-normal type) model as described in Durbin and Rocke
(2001,2002).
Then for each gene the mean of replicates for that gene is computed.
Optionally, if plot.it
is TRUE
then the mean
is plotted against its standard deviation (over
replicates).
Then the intensities are sorted according to increasing replicate mean.
Optionally, if plot.it
is TRUE
then a plot of the
intensities values as a vector (sorted according to increasing
replicate mean) is plotted in black, and then the true mean plotted
in colour 2 (on my screen this is red) and the computed replicate
mean plotted in green.
The DDHF transform of the sorted intensities is computed.
Optionally, if plot.it
is TRUE
then a plot of the
transformed means versus the transformed standard deviations is plotted.
Followed by a time series plot of the transformed sorted intensities.
These can be studied to see how well DDHF has done the transformation.
Then two smoothing methods are applied the the DDHF transformed data.
One method is translation invariant, Haar wavelet universal thresholding.
The other method is the classical smoothing spline. If plot.it
is
TRUE
then these smoothed estimates are plotted in different
colours.
Then the mean estimated intensity for each gene is computed and this
is returned as the first column of a two-column matrix (ansm
).
The second column is the true underlying mean. The object
hftssq
contains a measure of error between the estimated
and true gene means.