simulateSpectrum {Cardinal} | R Documentation |
Simulate mass spectra or complete MS imaging experiments, including a possible baseline, spatial and spectral noise, mass drift, mass resolution, and multiplicative variation, etc.
A number of preset imaging designs are available for quick-and-dirty simulation of images.
These functions are designed for small proof-of-concept examples and testing, and may not scale well to simulating larger datasets.
simulateSpectrum(n = 1L, peaks = 50L, mz = rlnorm(peaks, 7, 0.3), intensity = rlnorm(peaks, 1, 0.9), from = 0.9 * min(mz), to = 1.1 * max(mz), by = 400, sdpeaks = sdpeakmult * log1p(intensity), sdpeakmult = 0.2, sdnoise = 0.1, sdmz = 10, resolution = 1000, fmax = 0.5, baseline = 0, decay = 10, units=c("ppm", "mz"), representation = c("profile", "centroid"), ...) simulateImage(pixelData, featureData, preset, from = 0.9 * min(mz), to = 1.1 * max(mz), by = 400, sdrun = 1, sdpixel = 1, spcorr = 0.3, sptype = "SAR", representation = c("profile", "centroid"), units=c("ppm", "mz"), as = c("MSImagingExperiment", "SparseImagingExperiment"), BPPARAM = bpparam(), ...) addShape(pixelData, center, size, shape=c("circle", "square"), name=shape) presetImageDef(preset = 1L, nruns = 1, npeaks = 30L, dim = c(20L, 20L), peakheight = 1, peakdiff = 1, sdsample = 0.2, jitter = TRUE, ...)
n |
The number of spectra to simulate. |
peaks, npeaks |
The number of peaks to simulate. Not used if |
mz |
The theoretical m/z values of the simulated peaks. |
intensity |
The mean intensities of the simulated peaks. |
from |
The minimum m/z value used for the mass range. |
to |
The maximum m/z value used for the mass range. |
by |
The step-size used for the observed m/z-values of the profile spectrum. |
sdpeaks |
The standard deviation(s) for the distributions of observed peak intensities on the log scale. |
sdpeakmult |
A multiplier used to calculate |
sdnoise |
The standard deviation of the random noise in the spectrum on the log scale. |
sdmz |
The standard deviation of the mass error in the observed m/z values of peaks, in units indicated by |
resolution |
The mass resolution as defined by |
fmax |
The proportion of the maximum peak height to use when defining the mass resolution. |
baseline |
The maximum intensity of the baseline. Note that |
decay |
A constant used to calculate the exponential decay of the baseline. Larger values mean the baseline decays more sharply at the lower mass range of the spectrum. |
units |
The units for |
representation |
Should a profile spectrum be returned or only the centroided peaks? |
BPPARAM |
An optional instance of |
pixelData |
A |
featureData |
A |
preset |
A number indicating a preset image definition to use. |
nruns |
The number of runs to simulate for each condition. |
sdrun |
A standard deviation giving the run-to-run variance. |
sdpixel |
A standard deviation giving the pixel-to-pixel variance. |
spcorr |
The spatial autocorrelation. Must be between 0 and 1, where |
sptype |
The type of spatial autocorrelation. |
as |
The class of object to be returned. |
... |
Additional arguments to pass to |
dim |
The dimensions of the preset image. |
peakheight |
Reference intensities used for peak heights by the preset. |
peakdiff |
A reference intensity difference used for the mean peak height difference between conditions, for presets that simulate multiple conditions. |
sdsample |
A standard deviation giving the amount of variation from the true peak heights for this simulated sample. |
jitter |
Should random noise be added to the location and size of the shapes? |
center |
The center of the shape. |
size |
The size of the shape (from the center). |
shape |
What type of shape to add. |
name |
The name of the added column. |
The simulateSpectrum()
and simulateImage()
functions are used to simulate mass spectra and MS imaging experiments. They provide a great deal of control over the parameters of the simulation, including all sources of variation.
For simulateImage()
, the user should provide the design of the simulated experiment via matching columns in pixelData
and featureData
, where each column corresponds to different morphological substructures or differing conditions. These design data frames are returned in the metadata()
of the returned object for later reference.
A number of presets are defined by presetImageDef()
, which returns only the pixelData
and featureData
necessary to define the experiment for simulateImage()
. These can be referenced for help in understanding how to define experiments for simulateImage()
.
The preset images are:
1: a centered circle
2: a topleft circle and a bottomright square
3: two corner squares and a centered circle
4: a centered circle with conditions A and B in different runs
5: a topleft circle and a bottomright square with conditions A and B in different runs
6: two corner squares and a centered circle; the circle has conditions A and B in different runs
7: matched pairs of circles with conditions A and B within the same runs; includes reference peaks
8: matched pairs of circles inside squares with conditions A and B within the same runs; includes reference peaks
9: a small sphere inside a larger sphere (3D)
The addShape()
function is provided for convenience when generating the pixelData
for simulateImage()
, as a simple way of adding morphological substructures using basic shapes such as squares and circles.
For simulateSpectrum
, a list
with elements:
mz
: a numeric vector of the observed m/z values
intensity
: a numeric vector or matrix of the intensities
For simulateImage
, a MSImagingExperiment
or a SparseImagingExperiment
.
For addShape
, a new PositionDataFrame
with a logical column added for the corresponding shape.
For presetImageDef
, a list with two elements: the pixelData
and featureData
to be used as input to simulateImage()
.
Kylie A. Bemis
simulateSpectrum
,
simulateImage
register(SerialParam()) set.seed(1) # generate a spectrum s <- simulateSpectrum(1) plot(intensity ~ mz, data=s, type="l") # generate a noisy low-resolution spectrum with a baseline s <- simulateSpectrum(1, baseline=2, sdnoise=0.3, resolution=100) plot(intensity ~ mz, data=s, type="l") # generate a high-resolution spectrum s <- simulateSpectrum(1, peaks=100, resolution=10000) plot(intensity ~ mz, data=s, type="l") # generate an image x <- simulateImage(preset=1, npeaks=10, dim=c(10,10)) m <- mz(metadata(x)$design$featureData) image(x, mz=m[5]) plot(x, coord=c(x=3, y=3))