CorrelatedObservation#
- class scarlet2.CorrelatedObservation(data, psf=None, wcs=None, channels=None, renderer=None, name='', correlation_function=None, mask=None, native_scale=None, n_eff=None)[source]#
Bases:
ObservationContent and definition of an observation with pixel correlations
The noise model is still assumed to be Gaussian, but with correlations between pixels. The implementation computes the goodness of fit in Fourier space from the noise power spectrum to avoid the expensive computation of/with an inverse banded matrix in configuration space.
Create an observation with a correlated-noise Gaussian likelihood.
The noise model is given by the measured correlation_function; the power spectrum used by the likelihood is derived from it. In most cases it is easier to use
from_observation()(for a coadd already on the desired grid) orfrom_resampling()(to resample the data here), orstack_observations()to combine several correlated observations that share a grid.- Parameters:
data (
array) – Observed data, 2D or 3D (C, H, W).psf – As for
Observation.wcs – As for
Observation.channels – As for
Observation.renderer – As for
Observation.name – As for
Observation.correlation_function (
dict) – 2D pixel correlation function keyed by integer (dy, dx) lag, as returned byscarlet2.measure.correlation_function(). Required.mask (
array, optional) – Boolean array marking invalid pixels. Defaults to all-valid.native_scale (
astropy.units.Quantity, optional) – Angular size of a native, pre-resampling pixel, e.g.0.17 * u.arcsec. Set this if the data sit on a grid finer than their true resolution (resampled by an upstream pipeline or byfrom_resampling()). Only1 / (native_scale / pixel_scale) ** 2of the pixels are then independent, and chi^2 is restricted to the spatial frequencies the native sampling supports. Requires wcs. Mutually exclusive with n_eff.n_eff (
sequenceofint, optional) – Per-channel effective degrees of freedom, one entry per channel. Overrides the count derived from mask and native_scale. Used bystack_observations()to carry through the (possibly heterogeneous) counts of the stacked observations.
- property N#
Effective number of degrees of freedom of the noise model
For an observation without resampling this is the number of unmasked pixels, as in
Observation.N. When the correlations come from resampling onto a finer grid, the noise covariance is rank-deficient: there are only as many independent noise values as unmasked pixels in the original observation. Using the pixel count there would dilutegoodness_of_fit()and the likelihood normalization by the surplus modes, which carry no information.This is the sum of the per-channel
n_eff.
- __call__()#
Evaluate the model
- check_set_renderer(frame)#
Check existence of
renderer, or set it by callingmatch()- Parameters:
frame (
Frame) – The frame to match- Return type:
None
- correlation_function: dict#
Measured data correlation function, from py:meth:scarlet2.measure.correlation_function
- data: Array#
Observed data
- eval_chi_square_in_box_and_border(scene, border_width=3)#
Evaluate the weighted mean (weighted by the inverse variance weights) of the squared residuals for each source. Chi square is also computed for the perimeter outside the box of with border_width.
- Parameters:
scene (
Scene) – Scene containing the sourcesborder_width (
int) – width of the border around the source box
- Return type:
Dictofsources indicesandtheir corresponding Dictofresiduals insideandoutside source box.
- classmethod from_observation(obs, patch_size=50, maxlength=12, native_scale=None)[source]#
Create a
CorrelatedObservationfrom data that already sit on the desired grid.Use this for a coadd (drizzled, warped, stacked, …) whose pixel grid matches the model frame. The method finds a patch of size patch_size with as few sources as possible, measures the 2D pixel correlation function in it (out to maxlength), and turns that into the noise power spectrum used by the likelihood. Truncating the correlation function at maxlength biases the power spectrum, severely so once the correlation length approaches maxlength.
To resample the data onto the model grid here instead, use
from_resampling().- Parameters:
obs (
Observation) – Observation containing the data and weight map. The weights are used only to locate bad and bright pixels for the correlation measurement.patch_size (
int) – Linear size (in pixels) of the source-free patch used to measure the correlation function. Clamped to the image size; if it reaches the image size the whole image is used.maxlength (
int) – Maximum lag (in pixels) of the 2D correlation function. Large enough to cover the extent of the correlations, and small compared to patch_size so every lag is averaged over many pixel pairs (the assertion enforces4 * maxlength <= patch_size).native_scale (
astropy.units.Quantity, optional) – Angular size of a native, pre-resampling pixel, e.g.0.17 * u.arcsec. Set this only if an upstream pipeline already resampled the data onto a grid finer than their true resolution: the pixels are then correlated and not all independent, so chi^2 is restricted to the spatial frequencies the native sampling supports. SeeCorrelatedObservation.__init__().
- Return type:
- classmethod from_resampling(obs, to_frame, lanczos_order=9, resample_psf=True, patch_size=50, maxlength=12)[source]#
Create a
CorrelatedObservationby resampling obs onto the grid of to_frame.The Lanczos resampling correlates the pixels, so the result needs a correlated-noise likelihood. This method measures the correlation function of the resampled data.
If the spatial grid of to_frame is finer than obs, the resampled noise is rank-deficient: chi^2 will automatically be restricted to the frequencies the native (obs) sampling supports.
- Parameters:
obs (
Observation) – Observation to resample. Its pixel scale is taken as the native resolution.to_frame (
Frame) – Frame describing the desired spatial sampling, typically the model frame.lanczos_order (
int) – Order of the Lanczos resampling kernel.resample_psf (
bool, optional) – Whether to resample obs.psf to to_frame. Set to False only if the PSF is already sampled at the resolution of to_frame.patch_size (
int) – Passed tofrom_observation().maxlength (
int) – Passed tofrom_observation().
- Return type:
- get(name=None)#
Get parameter(s) from this module
- Parameters:
name (
str, optional) – Name of parameter. If not set, returns all parameters.- Returns:
requested data arrays for parameters
- Return type:
dict
- goodness_of_fit(model)#
Evaluate the goodness of the model fit to the data
For a Gaussian noise model, the gof is defined as the averaged squared deviation of the model from the data, scaled by the variance of the data, aka mean chi squared \(\frac{1}{N}\sum_i=1^N w_i (m_i - d_i)^2\) with inverse variance weights \(w_i\).
Up to a normalization, the gof is identical to
log-likelihood.- Parameters:
model (
array) – The (pre-rendered) predicted data cube, typically from evaluatingScene- Return type:
float
- log_likelihood(model)#
The logarithm the likelihood of
datagiven model- Parameters:
model (
array) – The (pre-rendered) predicted data cube, typically from evaluatingScene- Return type:
float
- mask: Array#
Mask for invalid pixels
- match(frame, renderer=None)[source]#
Construct the mapping between frame (from the model) and this observation frame
See
Observation.match(). In addition, this method determines whether chi^2 can be evaluated on the padded grid of the renderer.
- n_eff: tuple#
Per-channel effective number of degrees of freedom of the noise model, see
NOne entry per channel. For a channel without resampling this is its number of unmasked pixels; for a channel resampled onto a finer grid it is the smaller count of independent noise values on the native grid.
Nis the sum over channels.
- name: str#
Name to describe the observation
- property parameters#
Parameters defined for this module
- Returns:
name: (node, param) mapping for all parameters
- Return type:
dict
- render(model, **kwargs)#
Render model in the frame of this observation
- replace(name, value)#
Replace member attribuge name with value
- Parameters:
name (
str) – Name of member to replacevalue (
any) – Value to replace member with
- Returns:
The modified module.
- Return type:
- replace_subtree(where, value, parameters=None)#
Replace the subtree where(self) with value and carry over all other parameters
Unlike
replace(), value does not need to have the same structure as the subtree it replaces, e.g. an array morphology can be replaced by aSersicMorphology. Parameters defined for this module remain attached to the same nodes, except those inside the replaced subtree, which are dropped.- Parameters:
where (
callable) – Function module -> node that selects the subtree to replace, as inequinox.tree_at()value (
any) – Replacement for the subtree. If this module has a frame, attributes of value given in sky coordinates are converted to pixel coordinates of that frame.parameters (
callable, optional) – Function module -> None that is called with the new module inside of aParameterscontext to define parameters for the replacement, e.g. by callingParameteron nodes of value.
- Returns:
The modified module.
- Return type:
Examples
>>> sersic = SersicMorphology(n=jnp.array(2.), size=jnp.array(2.), shape=morph.shape) >>> def sersic_parameters(morph): >>> Parameter(morph.n, name="n:0", prior=dist.Uniform(0.5, 6)) >>> Parameter(morph.size, name="size:0", constraint=constraints.positive, stepsize=0.1) >>> >>> scene_ = scene.replace_subtree(lambda sc: sc.sources[0].morphology, sersic, sersic_parameters)
See also
- set(values)#
Set parameter(s) from this module with values
- Parameters:
values (
dict[str,jnp.array]) – values to replace parameters with, identified by their name- Returns:
new module with parameter(s) replaced by values
- Return type:
- weights: Array#
Statistical weights (usually inverse variance) for
log_likelihood()