Reference

API reference

All public names in gspice.__all__ are documented below except the submodule gspice.legacy, which has its own page.

Models and fitting

Fit a model of the spectra.

gspice.Model

class gspice.Model(
    mean: ArrayLike,
    cov: Iterable[ArrayLike | LowRank],
    n_train: Iterable[int | np.integer[Any]],
    *,
    edges: Iterable[float | np.floating[Any] | np.integer[Any]] = (),
    center: float | np.floating[Any] | np.integer[Any] | None = None,
    noise: bool | np.bool_ = False,
    weighted: bool | np.bool_ = False,
    blocks: Blocks | None = None,
)

Fitted Gaussian model of spectra on a common pixel grid.

A spectrum in S/N bin $b$ has the working vector $x \sim N(0, \mathrm{cov}[b])$ in scaled units. Let $\tilde f$ be the flux of the spectrum with its bad pixels interpolated. Let $a$ be the square root of the mean inverse variance of its good pixels (see scale). Then $x = a (\tilde f - \mathrm{center}) - \mathrm{mean}[b]$ for a constant center, and $x = a \tilde f - \mathrm{mean}[b]$ otherwise.

A noise-aware model (noise=True) has one bin and works in flux units. Spectrum $s$ has

x = \tilde f - \mathrm{mean}[0] \sim N(0, C_s), \qquad C_s = \mathrm{cov}[0] + \mathrm{diag}(\nu_s),

where cov[0] is the covariance of the noise-free flux. $\nu_{si}$ is $1/\mathrm{ivar}_{si}$ at the good pixels and $1/a^2$ at the bad pixels of the spectrum. The constructor checks only the diagonal of cov[0]. Prediction rejects a $C_s$ that is not positive definite where it factors $C_s$ . The prior factors only $C_s$ at the target pixels.

The constructor takes the attributes below as array_like values, iterables and numpy scalars. It stores them as the types given here, with read-only float64 copies of the arrays and new LowRank covariances. The attributes after n_train are keyword-only. The constructor does not check that a covariance is positive definite, which costs $\mathrm{npix}^3$ . Two models are equal only if they are the same object.

A model computes the precision matrix of a bin at its first use and keeps it. It also keeps the guard weights of the latest nguard of each bin (see pixelwise). Each matrix takes $\mathrm{npix}^2 \times 8$ bytes, 46 MB at npix 2400, also for a LowRank covariance. Pickle, copy and deepcopy start without them.

Attributes:

mean (ndarray, shape (nbin, npix)) – The vector mean[b] of each of the nbin S/N bins (scaled; flux for a noise-aware model), subtracted to give $x$ , over the npix pixels. Finite; nbin and npix are 1 or more.
cov (tuple of ndarray, each shape (npix, npix), or LowRank) – The covariance of $x$ in each bin (scaled²), or of the noise-free flux of a noise-aware model (flux²). Either nbin finite, exactly symmetric matrices with a positive diagonal, or nbin LowRank covariances of npix pixels with positive diag. All bins hold one kind. A noise-aware model allows a diagonal ≥ 0, or a LowRank with diag ≥ 0. A positive semi-definite cov[0] makes every $C_s$ positive definite.
n_train (tuple of int) – The number $N_s$ of training spectra in each bin (nbin integers ≥ 1).
edges (tuple of float) – The nbin − 1 interior S/N edges between the bins: finite and strictly increasing. Spectrum $s$ with the signal-to-noise ratio $\mathrm{snr}_s$ (see snr) is in the bin np.searchsorted(edges, snr_s, side="right").
center (float or None) – None for a mean-centered model, or the constant center $c$ (flux). None for a noise-aware model.
noise (bool) – True for a noise-aware model, which has one bin and center None.
weighted (bool) – True if the model was fitted with positive weights that are not all equal in a bin.
blocks (Blocks or None) – The blocks and thresholds of fit_blocks, with a row of thresholds for each bin and every block within the npix pixels. None for a model without them.

Raises:

TypeError – If a field has a type or dtype that is not allowed above. The message starts with the name of the field.
ValueError – If a field has a shape or a value that is not allowed above. The message starts with the name of the field.

gspice.LowRank

class gspice.LowRank(vectors: ArrayLike, values: ArrayLike, diag: ArrayLike)

Low-rank plus diagonal covariance.

C = V \mathrm{diag}(\mathrm{values}) V^\mathsf{T} + \mathrm{diag}(\mathrm{diag})

fit with rank=K returns a Model with a LowRank covariance. A model computes its precision matrix by the Woodbury identity, without forming $C$ . It rejects the covariance where $\kappa_W > 10^{13}$ , with

\kappa_W = \frac{\max \tilde d + \lambda_{\max}\bigl( \mathrm{diag}(\mathrm{values})^{1/2} \tilde V^\mathsf{T} \tilde V \, \mathrm{diag}(\mathrm{values})^{1/2}\bigr)}{\min \tilde d}.

Here $c_i = \mathrm{diag}_i + \sum_\ell \mathrm{values}_\ell V_{i\ell}^2$ are the variances of $C$ , $\tilde d = \mathrm{diag} / c$ and $\tilde V = \mathrm{diag}(c)^{-1/2} V$ . $\kappa_W$ bounds both the condition number of $C$ scaled to unit diagonal and $\max_i c_i / \mathrm{diag}_i$ . The constructor takes array_like values and stores read-only float64 copies. Two LowRank are equal only if they are the same object.

Attributes:

vectors (ndarray, shape (npix, K)) – The $K \ge 1$ columns of $V$ , for npix ≥ 1 pixels: finite. fit returns orthonormal columns (unit-free), in decreasing order of values.
values (ndarray, shape (K,)) – The variance of each column (scaled²; flux² for a noise-aware model): finite and positive.
diag (ndarray, shape (npix,)) – The diagonal term (scaled²; flux² for a noise-aware model): finite and 0 or more. Model requires positive values, except for a noise-aware model.

Raises:

TypeError – If a field has a dtype that is not integer or floating. The message starts with the name of the field.
ValueError – If a field has a shape or a value that is not allowed above. The message starts with the name of the field.

gspice.Blocks

class gspice.Blocks(
    starts: ArrayLike,
    width: int | np.integer[Any],
    nguard: int | np.integer[Any],
    thresholds: ArrayLike,
)

Geometry of the block scan and threshold of each block, from fit_blocks.

Block $j$ holds the pixels starts[j] to starts[j] + width - 1. Its statistic $\chi^2_A$ of fit_blocks compares the block with its prediction from the pixels farther than nguard from it. thresholds[b, j] is the threshold of $\chi^2_A$ for the spectra of S/N bin $b$ . The constructor takes array_like values, of an integer dtype for starts, and numpy scalars. It stores them as the types given here, with read-only copies of the arrays. Two Blocks are equal only if they are the same object.

Attributes:

starts (ndarray of int64, shape (nblock,)) – The first pixel of each of the nblock ≥ 1 blocks: 0 or more, and strictly increasing.
width (int) – The number of pixels per block, $|A|$ ; 1 or more.
nguard (int) – The half-width of the guard window around each block (pixels); from 0 to $2^{63} - 1$ , so that save can write it.
thresholds (ndarray, shape (nbin, nblock)) – The threshold of $\chi^2_A$ of each bin and block (unit-free): finite and positive.

Raises:

TypeError – If a field has a type or dtype that is not allowed above. The message starts with the name of the field.
ValueError – If a field has a shape or a value that is not allowed above. The message starts with the name of the field.

gspice.fit

gspice.fit(
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
    *,
    nsigma: Sequence[float | np.floating[Any] | np.integer[Any]] | NDArray[np.integer[Any] | np.floating[Any]] = (20.0, 8.0, 6.0),
    maxbadpix: int | np.integer[Any] = 64,
    nguard: int | np.integer[Any] = 20,
    cumulative: bool | np.bool_ = True,
    center: float | np.floating[Any] | np.integer[Any] | None = None,
    weights: ArrayLike | None = None,
    rank: int | np.integer[Any] | None = None,
    noise: bool | np.bool_ = False,
    snr_bins: int | np.integer[Any] | Iterable[float | np.floating[Any] | np.integer[Any]] | None = None,
    progress: Callable[[str, int, int], object] | None = None,
) -> tuple[Model, NDArray[np.bool_]]

Fit a model to training spectra with iterative outlier clipping (paper Alg. 4).

The bad pixels are (mask != 0) | (ivar == 0), and fit never reads flux or ivar there. fit uses a spectrum unless it has more than maxbadpix bad pixels, has the weight 0, or a pass mask or the final mask covers it completely. $N_s$ counts the spectra of a fit. A fit on a mask interpolates the pixels in it, and scales each spectrum by $a_s$ over the pixels outside it (see scale). It then estimates

m = \frac{1}{V_1} \sum_s w_s x_s, \qquad S = \frac{\sum_s w_s (x_s - m)(x_s - m)^\mathsf{T}}{V_1 - V_2 / V_1},

where $x_s = a_s \tilde f_s$ , $V_1 = \sum_s w_s$ and $V_2 = \sum_s w_s^2$ . With unit weights, $m$ is the column mean and the divisor is $N_s - 1$ . With a constant center $c$ , $m = 0$ , $x_s = a_s (\tilde f_s - c)$ and $S = \sum_s w_s x_s x_s^\mathsf{T} / V_1$ . With rank=K, the covariance is the LowRank $C = V \mathrm{diag}(\lambda) V^\mathsf{T} + \mathrm{diag}(d)$ . Here $(V, \lambda)$ are the $K$ leading eigenpairs of $S$ , and $d = \mathrm{diag}(S) - \sum_\ell \lambda_\ell V_\ell^2$ , so that $C$ and $S$ have the same diagonal.

Pass $n$ of the $N$ = len(nsigma) passes fits on its mask $B_n$ with unit weights and the rank, and flags $F_n$ : each pixel, interpolated ones included, whose in-sample $|z|$ of pixelwise with nguard exceeds nsigma[n - 1], and its two neighbors. $B_1$ holds the bad pixels.

With cumulative=True (paper Alg. 4), $B_{n+1} = B_n \cup F_n$ , and the final mask is $B_{N+1}$ .
With cumulative=False (2022), $B_{n+1} = B_1 \cup F_n$ , and the final mask is $B_1 \cup F_{N-1} \cup F_N$ .
nsigma=() fits on the mask of the bad pixels, without clipping.

The final fit uses the final mask and the weights $w_s$ . Without rank, fit logs a warning on the "gspice" logger for each pass whose threshold exceeds the in-sample ceiling of $|z|$ of an outlier in a spectrum of typical leverage.

With noise=True, the passes clip as above. The final fit then estimates the covariance of the noise-free flux, in flux units. Let $u_{si}$ be 1 at the pixels outside the final mask and 0 elsewhere, and $v_s = a_s^2 w_s$ . Then

\begin{aligned} W_i &= \sum_s v_s u_{si}, \qquad \mu_i = \frac{1}{W_i} \sum_s v_s u_{si} f_{si}, \qquad A_{si} = v_s^{1/2} u_{si} (f_{si} - \mu_i), \\ D^{(1)}_{ij} &= \sum_s v_s u_{si} u_{sj}, \qquad D^{(2)}_{ij} = \sum_s v_s^2 u_{si} u_{sj}, \qquad N_i = \sum_s \frac{v_s u_{si} (1 - v_s / W_i)}{\mathrm{ivar}_{si}}, \\ \Sigma_{ij} &= \frac{(A^\mathsf{T} A)_{ij} - \delta_{ij} N_i} {D^{(1)}_{ij} - D^{(2)}_{ij} \bigl(1/W_i + 1/W_j - D^{(1)}_{ij} / (W_i W_j)\bigr)}. \end{aligned}

$\Sigma$ is unbiased element by element under these conditions:

The spectra are independent, with a common mean and signal covariance.
The noise has zero mean and the variance $1/\mathrm{ivar}$ . It is uncorrelated between pixels and with the signal.
The masks, weights and ivar do not depend on the flux. The outlier flags of the final mask break this condition.

The model holds $\mathrm{mean}[0] = \mu$ and $\mathrm{cov}[0] = Q \max(\Lambda, 0) Q^\mathsf{T}$ , with $\Sigma = Q \Lambda Q^\mathsf{T}$ its eigendecomposition. This cov[0] is biased upward. With rank=K, cov[0] is the LowRank $C = V \mathrm{diag}(\lambda) V^\mathsf{T} + \mathrm{diag}(d)$ . Here $(V, \lambda)$ are the $K$ leading eigenpairs of $\Sigma$ , and $d_i = \sum_{\ell > K} \max(\lambda_\ell, 0) Q_{i\ell}^2$ . So $C$ and $Q \max(\Lambda, 0) Q^\mathsf{T}$ have the same diagonal.

With snr_bins, spectrum $s$ goes to the S/N bin np.searchsorted(edges, snr_s, side="right"), with $\mathrm{snr}_s$ from snr on flux, ivar and mask. Each bin runs the passes and the final fit above on its own spectra. Bin $b$ of the model thus equals a fit on the spectra of bin $b$ alone, with their own clipping, $m$ , $S$ and $N_s$ . Each bin clips with its own $N_s$ . So the in-sample ceiling of $|z|$ falls in small bins, and they flag fewer outliers.

With cumulative=False, fit equals gspice.legacy.covar_iter_mask: the covariance to rounding, and the mask exactly. This holds under these conditions:

nguard is 20, center, weights, rank and snr_bins are None, and noise is False.
flux is float64, and ivar == 0 only at pixels in mask.
nsigma is not empty.
No $|z|$ lies within the error of the 2022 Schur updates of a threshold.
No pass mask or final mask covers completely a spectrum with at most maxbadpix bad pixels.
fit raises no error.

Parameters:

flux, ivar, mask (array_like, shape (nspec, npix)) – As in scale, except that fit drops a spectrum without a good pixel.
nsigma (sequence of float or ndarray, shape (N,), optional) – The threshold on $|z|$ of each clipping pass, in order; finite and positive.
maxbadpix (int, optional) – The largest number of bad pixels of a used spectrum; 0 or more.
nguard (int, optional) – The half-width of the guard window of z (pixels); 0 or more.
cumulative (bool, optional) – Whether each pass keeps the flags of every earlier pass (paper Alg. 4), or only those of the pass before it (2022).
center (float, optional) – The constant center $c$ (flux); finite. None: a mean-centered model.
weights (array_like, shape (nspec,), optional) – The weights $w_s$ of the final fit, of an integer or floating dtype; finite and 0 or more. None: equal weights.
rank (int, optional) – The rank $K$ $K$ of a low-rank model, used by every pass; from 1 to npix −
1. None: a dense model.
noise (bool, optional) – If True, fit a noise-aware model (see Model). It needs center and snr_bins None.
snr_bins (int or sequence of float, optional) – The S/N bins. Either their number $n$ , 1 or more, with the interior edges at the $1/n, \dots, (n - 1)/n$ quantiles (np.quantile, linear) of $\mathrm{snr}_s$ of the spectra that fit uses. These spectra have at most maxbadpix bad pixels, a positive weight and a good pixel. Or the interior edges, finite and strictly increasing. None: one bin.
progress (callable, optional) – progress(stage, done, total), called after each chunk of spectra of each pass and of the final fit. Its arguments are the stage "fit", the number of spectra of that fit done, and their total $N_s$ . With S/N bins, each fit of a bin counts the spectra of that bin.

Returns:

model (Model) – The model, with the edges, the $N_s$ of each bin in n_train, and the center. weighted is True if the weights of the spectra of a bin in the final fit are not all equal. The covariances are LowRank with rank=K. The model holds its precision matrices (see Model), except a noise-aware model, whose prediction forms one per spectrum. For a noise-aware model, weighted is True if the weights $w_s$ are not all equal, whatever $a_s$ .
mask (ndarray of bool, shape (nspec, npix)) – The final mask, and True at every pixel of a spectrum that fit does not use.

Raises:

TypeError – If an argument has a type or a dtype that is not allowed above.
ValueError – If an argument has a shape or a value that is not allowed above. If a fit, of a bin with snr_bins, lacks the spectra that it needs:
- For a dense model, $N_s - 1 \ge \mathrm{npix}$ if it is mean-centered, and $N_s \ge \mathrm{npix}$ for a constant center. Also a good pixel in one of them at every pixel.
- With rank=K, $N_s - 1 \ge K + 1$ if it is mean-centered, and $N_s \ge K + 1$ for a constant center, for $S$ of rank $K + 1$ or more.
With rank=K, also if $\lambda_K \le 0$ or a $d_i \le 0$ . If noise is True with a center or snr_bins. With noise=True, the passes need the spectra above. The final fit raises instead if fewer than two spectra are good at a pixel, or if a denominator of $\Sigma$ is 0 or negative. With rank=K, it also raises if $\lambda_K \le 0$ for $\Sigma$ . With two spectra at each pixel, a denominator is 0 or negative where no spectrum is good at both pixels of a pair. Rounding near 0 can also make it 0 or negative. The message names the pass or the final fit, the bin and its edges, and the counts. It also names the first such pixels or pairs of pixels, or $\lambda_K$ . fit checks the counts of every bin before any bin clips.
numpy.linalg.LinAlgError – If a covariance has a variance of 0, or fails the check of its precision matrix. The check fails if the matrix is not positive definite. It also fails if the reciprocal condition number of the covariance, scaled to unit diagonal, is below $10^{-13}$ . This can happen for $N_s - 1$ close to npix. With rank=K, the check fails if $\kappa_W$ of LowRank is above $10^{13}$ . Clipping checks $P$ only if nguard < npix - 1. The message names the pass or the final fit, the bin and $N_s$ . It also names the pixel of variance 0, the reciprocal condition number below $10^{-13}$ , or $\kappa_W$ . With noise=True, only for the clipping covariances.

Notes

Cost:

Per pass and bin, with the $N_s$ of the bin: about $N_s \, \mathrm{npix}^2$ flops for the covariance, and $\mathrm{npix}^3$ for its precision matrix. The guard weights cost $\mathrm{npix} (w^3/3 + 2 \, \mathrm{npix} \, w)$ , and z costs $2 N_s \, \mathrm{npix}^2$ . Here $w = \min(2 \, \mathrm{nguard} + 1, \mathrm{npix})$ .
The final fit: about $N_s \, \mathrm{npix}^2 + \mathrm{npix}^3$ flops per bin.
With nbin above 1, the S/N sorts the npix values of each spectrum once.
With rank=K, each fit adds the $K$ leading eigenpairs of $S$ . Where $K \le \mathrm{npix}/80$ , Lanczos costs about $2 \, \mathrm{npix}^2$ flops per product with $S$ . Otherwise they cost about $(4/3) \, \mathrm{npix}^3$ flops, more than the $\mathrm{npix}^3$ of a dense $P$ . $P$ then costs about $\mathrm{npix}^2 K$ flops.
Each fit holds the working vectors of its spectra, $N_s \, \mathrm{npix} \times 8$ bytes.
With noise=True, the final fit costs about $3 N_s \, \mathrm{npix}^2$ flops for $A^\mathsf{T} A$ , $D^{(1)}$ and $D^{(2)}$ . It also costs a full eigendecomposition of an npix × npix matrix, also with rank=K. It holds a few arrays of one chunk of max(1, 2**24 // npix) spectra instead of the working vectors.

gspice.fit_blocks

gspice.fit_blocks(
    model: Model,
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
    *,
    width: int | np.integer[Any] = 150,
    stride: int | np.integer[Any] = 75,
    nguard: int | np.integer[Any] = 20,
    percentile: float | np.floating[Any] | np.integer[Any] = 99.9,
    progress: Callable[[str, int, int], object] | None = None,
) -> Model

Return the model with block thresholds from clean spectra.

Block $j$ holds the width pixels from starts[j]. The starts are 0, stride, 2 stride and so on, up to npix − width. If the last of them is smaller than npix − width, the start npix − width follows. The defaults give 31 blocks at npix 2400. fit_blocks first replaces the bad pixels by their conditional mean, as infill does. It does not detect outliers by z.

Take a spectrum in S/N bin $b$ , with $x$ the working vector of its infilled flux (see Model) and $P$ the inverse of cov[b]. Let $A$ be a block, $G$ the pixels within nguard of $A$ , clipped at the ends, and $O$ the pixels outside $G$ . Then

\chi^2_A = \frac{r_A^\mathsf{T} C_{AA}^{-1} r_A}{|A|}, \qquad r_A = x_A - \mu_A,

where $\mu_A$ and $C_{AA} = (P_{GG}^{-1})_{AA}$ are the conditional mean and covariance of $x_A$ given $x_O$ . The guard pixels $G \setminus A$ are thus marginalized. Where $G$ covers all pixels, $\mu_A$ and $C_{AA}$ are those of the prior, 0 and cov[b] at $A$ , computed without $P$ .

thresholds[b, j] is np.percentile (linear) of $\chi^2_A$ of block $j$ , at $p$ = percentile, over the $m$ spectra of bin $b$ . A block without an artifact then exceeds its threshold with the probability

\frac{m (1 - p/100) + p/100}{m + 1}

if $(m - 1) p/100$ is an integer. Otherwise the probability is at most this value, for a density that falls beyond the threshold. For $p = 99.9$ , it is 0.20% at $m = 1000$ and 0.12% at $m = 5000$ . repair also infills the pixels that z flags, which lowers this rate. Use clean held-out spectra.

For a noise-aware model, $x$ is in flux units. The covariance $C_s$ of each spectrum (see Model) then replaces cov[b] in its infill and its blocks.

Parameters:

model (Model) – The model.
flux, ivar, mask (array_like, shape (nspec, npix)) – As in scale, with the npix of model. Clean spectra, with at least $m$ spectra in each bin such that $m (100 - p) \ge 100$ . The check reads $p$ as the decimal Fraction(str(p)): 1000 spectra at $p = 99.9$ .
width (int, optional) – The number of pixels per block, $|A|$ ; from 1 to npix.
stride (int, optional) – The distance between the starts of consecutive blocks (pixels); 1 or more.
nguard (int, optional) – The half-width of the guard window around each block (pixels); from 0 to $2^{63} - 1$ , as for Blocks.
percentile (float, optional) – The percentile $p$ of the thresholds; $0 < p < 100$ .
progress (callable, optional) – progress(stage, done, total), called after each chunk of spectra, with the stage "fit_blocks", the number of spectra done and their total. For a noise-aware model, a chunk is one spectrum.

Returns:

Model – A new model, equal to model except for blocks: a Blocks that holds the starts, width, nguard and thresholds.

Raises:

TypeError – If model is not a Model, or width, stride or nguard is not an integer. If percentile is not a real number, or progress is not callable or None. If flux, ivar or mask has a dtype that scale does not allow.
ValueError – If flux, ivar or mask has a shape or a value that scale does not allow. If flux has another npix than model, or a spectrum has no good pixel. If width, stride, nguard or percentile is outside its range. If a bin has fewer spectra than $p$ needs; the message names the bin and its count. If a threshold is not finite or is 0, as for spectra equal to their conditional mean.
numpy.linalg.LinAlgError – If the covariance of a bin fails the check of its precision matrix, as in pixelwise, and the $G$ of a block does not cover all pixels. If $P_{GG}$ of a block is not positive definite. If cov[b] at $T = A$ fails that check, for a block whose $G$ covers all pixels. Also as in infill, for the bad pixels. The message names the bin, and the block or the row of a pattern. For a noise-aware model, as in infill, with $C_s$ in place of cov[b].

Notes

Cost:

$\mathrm{npix}^3$ flops for $P$ at the first use of a bin, then $|G|^3/3$ flops per bin and block.
About $2\,\mathrm{npix}^2 + \sum_A |G|^2$ flops per spectrum, plus the infill of the bad pixels, of the mean only.
The factors of all bins take $\mathrm{nbin} \sum_A |G|^2 \times 8$ bytes during the scan: 8.8 MB per bin at npix 2400 with the defaults, and 649 MB at stride=1.
Beyond the inputs, and their float64 copies if they are long double, fit_blocks holds a bool array of the input shape. It also holds the $\chi^2_A$ of the spectra of one bin (nblock float64 per spectrum), and a few float64 arrays of one chunk of max(1, 2**24 // npix) spectra.
The new model shares the arrays of model, its precision matrices and its guard weights (see Model).
Unless every block covers all pixels with its guard, a noise-aware model forms $P_s$ and the factors of the blocks for each spectrum. $P_s$ takes $\mathrm{npix}^3$ flops for a dense cov[0], or $\mathrm{npix}^2 K$ for a LowRank of rank $K$ . Otherwise the infill of the bad pixels forms no $P_s$ , and costs as in infill.

Prediction

Predict pixels from other pixels.

gspice.Prediction

class gspice.Prediction(
    flux: NDArray[np.float64],
    var: NDArray[np.float64],
    z: NDArray[np.float64],
)

The prediction of each pixel of a set of spectra, from pixelwise.

The arrays that pixelwise returns are read-only.

Attributes:

flux (ndarray, shape (nspec, npix)) – The prediction of each pixel (flux).
var (ndarray, shape (nspec, npix)) – The conditional variance of each pixel (flux²), which includes the noise. With inflate=True, it is multiplied by $f_s(N_s, N_k, 1)$ .
z (ndarray, shape (nspec, npix)) – $(f - \mathrm{flux}) / \sqrt{\mathrm{var}}$ at good pixels, with $f$ the input flux, and NaN at bad pixels (unit-free).

gspice.Repaired

class gspice.Repaired(
    flux: NDArray[np.float64],
    var: NDArray[np.float64],
    replaced: NDArray[np.bool_],
    ignored: NDArray[np.bool_],
    chi2: NDArray[np.float64],
    dof: NDArray[np.int64],
    fs: NDArray[np.float64],
    draws: NDArray[np.float64] | None,
    blocks: NDArray[np.bool_] | None = None,
    _covariances: _Covariances | None = None,
)

Spectra whose replaced pixels hold their conditional mean, from infill.

repair also returns one. The arrays that infill and repair return are read-only. With covariance=True, covariance(s) gives the conditional covariance of the replaced pixels of spectrum s.

Attributes:

flux (ndarray, shape (nspec, npix)) – The conditional mean at the replaced pixels, and the input flux elsewhere (flux).
var (ndarray, shape (nspec, npix)) – The conditional variance at the replaced pixels, which includes the noise. With inflate=True, it is multiplied by $f_s(N_s, N_k, 1)$ . 1/ivar at the other good pixels, and inf at the other bad pixels (flux²).
replaced (ndarray of bool, shape (nspec, npix)) – The replaced pixels $M$ .
ignored (ndarray of bool, shape (nspec, npix)) – The pixels that are neither conditioned on nor replaced, $G \setminus M$ .
chi2 (ndarray, shape (nspec,)) – $r_T^\mathsf{T} C_{TT}^{-1} r_T$ over the replaced pixels $T$ that are good in the input (unit-free). $r_T$ is their input values minus the conditional mean. $C_{TT}$ is their conditional covariance, which inflate=True multiplies by $f_s(N_s, N_k, 1)$ . 0 where dof is 0.
dof (ndarray of int64, shape (nspec,)) – The number $|T|$ of replaced pixels that are good in the input.
fs (ndarray, shape (nspec,)) – $f_s(N_s, N_k, \mathrm{dof})$ of finite_sample_factor (unit-free), divided by $f_s(N_s, N_k, 1)$ with inflate=True. $N_s$ is the n_train of the bin of the spectrum, and $N_k = \mathrm{npix} - |G|$ the number of its conditioning pixels. NaN where dof is 0 or $N_s - N_k - \mathrm{dof} - 2 \le 0$ . 1 for a model with a constant center, weights, a LowRank covariance or noise, for which no factor is derived. Take the other models, and replace and ignore masks that do not depend on the data. Then chi2 / (dof * fs) has the mean 1 for held-out Gaussian spectra, with or without inflate.
draws (ndarray, shape (nspec, npix), or None) – A draw from the conditional distribution at the replaced pixels, with the covariance of covariance, and the input flux elsewhere (flux). None if the call ran without rng.
blocks (ndarray of bool, shape (nspec, nblock), or None) – The blocks of model.blocks that the block pass of repair flagged in each spectrum. None for a result of infill, or of repair with blockwise=False.

covariance

covariance(s: int | np.integer[Any]) -> tuple[NDArray[np.int64], NDArray[np.float64]]

Return the conditional covariance of the replaced pixels of spectrum s.

Parameters:

s (int) – The row of the spectrum; from 0 to nspec − 1.

Returns:

pixels (ndarray of int64, shape (m,)) – The replaced pixels of spectrum s, in increasing order.
cov (ndarray, shape (m, m)) – Their conditional covariance (flux²), multiplied by $f_s(N_s, N_k, 1)$ with inflate=True. Its diagonal is var[s, pixels].

Raises:

TypeError – If s is not an integer.
ValueError – If the call ran without covariance=True, or s is outside 0 to nspec − 1.

gspice.pixelwise

gspice.pixelwise(
    model: Model,
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
    *,
    nguard: int | np.integer[Any] = 20,
    inflate: bool | np.bool_ = False,
    progress: Callable[[str, int, int], object] | None = None,
) -> Prediction

Predict each pixel from all pixels farther than nguard from it (paper Alg. 3).

Take pixel $i$ of a spectrum in S/N bin $b$ , with $x$ its working vector (see Model) and $P$ the inverse of cov[b]. Let $g = \{j : |j - i| \le \mathrm{nguard}\}$ , and let $O$ hold the pixels outside $g$ . Then

x_i \mid x_O \sim N\bigl(-(P_{gg}^{-1} P_{gO} x_O)_i, \ (P_{gg}^{-1})_{ii}\bigr).

Prediction.flux and Prediction.var are this mean and variance in flux units, as unscale returns them. Bad pixels enter as interpolated values. Their z is NaN, and var omits the error of those values. Next to a run of 60 bad pixels at the correlation length 20, the measured $\langle z^2 \rangle$ was 0.65 to 2.45, depending on the noise level.

A bad pixel farther than nguard from $i$ takes its value from the nearest good pixels. So a good pixel within nguard of $i$ can enter the prediction of $i$ . Next to a run of bad pixels longer than nguard, $i$ itself can also enter it. nguard >= npix - 1 gives the prior.

For a noise-aware model, $x$ is in flux units, and $P$ is the inverse of the covariance $C_s$ of the spectrum (see Model). pixelwise forms $P$ for each spectrum. Where $2 \, \mathrm{nguard} + 1 \le \mathrm{npix}$ , it computes the mean as $x_i - (P_{gg}^{-1} (P x)_g)_i$ . This is equal in exact arithmetic, but the values in $g$ cancel only to rounding.

Parameters:

model (Model) – The model.
flux, ivar, mask (array_like, shape (nspec, npix)) – As in scale, with the npix of model.
nguard (int, optional) – The half-width of the guard window $g$ (pixels); 0 or more.
inflate (bool, optional) – If True, multiply var by $f_s(N_s, N_k, 1)$ of finite_sample_factor. $N_s$ is the n_train of the bin of the spectrum. $N_k = \mathrm{npix} - |g|$ is the number of conditioning pixels of each pixel, which is larger near the ends, where $g$ is clipped. $z^2$ then has the mean 1 for held-out Gaussian spectra. It needs a dense, mean-centered, unweighted scaled model.
progress (callable, optional) – progress(stage, done, total), called after each chunk of spectra, with the stage "pixelwise", the number of spectra done and their total. For a noise-aware model, a chunk is one spectrum.

Returns:

Prediction – The prediction, in flux units.

Raises:

TypeError – If model is not a Model, nguard is not an integer, inflate is not a bool, or progress is not callable or None. If flux, ivar or mask has a dtype that scale does not allow.
ValueError – If flux, ivar or mask has a shape or a value that scale does not allow. If flux has another npix than model, nguard is negative, or a spectrum has no good pixel. If inflate is True and model has a constant center, weights, a LowRank covariance or noise. If inflate is True and $N_s - N_k - 3 \le 0$ at a pixel of a bin with spectra, where $f_s(N_s, N_k, 1)$ is not defined; the message names $N_s$ and $N_k$ .
numpy.linalg.LinAlgError – If nguard < npix - 1 and the covariance of a bin with spectra fails the check of its precision matrix. The check fails if the matrix is not positive definite, or if the reciprocal condition number of the covariance, scaled to unit diagonal, is below $10^{-13}$ . It also fails if $\kappa_W$ of a LowRank covariance is above $10^{13}$ . The message names the bin, or for a noise-aware model the row whose $C_s$ fails.

Notes

Cost:

$\mathrm{npix}^3$ flops for $P$ at the first use of a bin, or $\mathrm{npix}^2 K$ for a LowRank of rank $K$ .
About $\mathrm{npix} (w^3/3 + 2 \, \mathrm{npix} \, w)$ flops for the guard weights, at the first use of a bin and after each change of nguard, with $w = \min(2 \, \mathrm{nguard} + 1, \mathrm{npix})$ .
Then $2 \, \mathrm{npix}^2$ flops per spectrum.
The model keeps $P$ and the guard weights (see Model).
A noise-aware model costs $P$ per spectrum, and keeps none. Where $2 \, \mathrm{nguard} + 1 \le \mathrm{npix}$ , it costs about $\mathrm{npix} (w^3/3 + 2 w^2)$ flops for its guard windows and $2 \, \mathrm{npix} (\mathrm{npix} + w)$ for the mean, instead of the guard weights.

gspice.infill

gspice.infill(
    model: Model,
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
    *,
    replace: ArrayLike | None = None,
    ignore: ArrayLike | None = None,
    inflate: bool | np.bool_ = False,
    covariance: bool | np.bool_ = False,
    rng: np.random.Generator | None = None,
    progress: Callable[[str, int, int], object] | None = None,
) -> Repaired

Replace pixels by their Gaussian conditional mean given the conditioning pixels.

Take a spectrum in S/N bin $b$ , with $x$ its working vector (see Model) and $P$ the inverse of cov[b]. Let $G$ be the pixels of replace | ignore | bad, with the bad pixels (mask != 0) | (ivar == 0). The other pixels $O$ are the conditioning pixels. Then

x_G \mid x_O \sim N(-P_{GG}^{-1} P_{GO} x_O, \ P_{GG}^{-1}).

Repaired.flux is this mean at the replaced pixels $M$ = replace, in flux units as unscale returns it, and the input flux elsewhere. The pixels of $G$ outside $M$ keep their input value, and the result reports them as ignored. A pixel in both replace and ignore is replaced. Values in $G$ never enter the mean, so an outlier there cannot change it. chi2 compares the input values at $T = M \setminus \mathrm{bad}$ , the replaced pixels that are good in the input, with their conditional distribution.

Where $G$ covers all pixels, the result is the prior $N(0, \mathrm{cov}[b])$ , which infill computes from cov[b] without $P$ . For a noise-aware model, $x$ is in flux units, and the $C_s$ of the spectrum (see Model) replaces cov[b]. Each spectrum then has its own $P$ and pattern.

With rng, Repaired.draws holds the draw $\mu_M + F z_M$ of $x_M$ , in flux units like the mean. Here $\mu_M$ is the mean at $M$ , $F F^\mathsf{T} = (P_{GG}^{-1})_{MM}$ , and $z_M$ holds the normals of the spectrum at $M$ .

Spectra of one S/N bin with equal $G$ , $M$ and $T$ form a pattern. The spectra of a pattern in one chunk of max(1, 2**24 // npix) rows share one factorization of $P_{GG}$ . To infill with another model the pixels that a result r replaced and ignored, pass replace=r.replaced and ignore=r.ignored.

Parameters:

model (Model) – The model.
flux, ivar, mask (array_like, shape (nspec, npix)) – As in scale, with the npix of model.
replace (array_like of bool, shape (nspec, npix) or (npix,), optional) – The pixels to replace, $M$ . An (npix,) mask applies to every spectrum, for example to predict a range of pixels. None: the bad pixels.
ignore (array_like of bool, shape (nspec, npix) or (npix,), optional) – Pixels that are neither conditioned on nor replaced, unless they are in replace, such as a guard around a range of pixels. None: no pixels.
inflate (bool, optional) – If True, multiply var and the covariance at the replaced pixels by $f_s(N_s, N_k, 1)$ of finite_sample_factor. $N_s$ is the n_train of the bin of the spectrum, and $N_k = \mathrm{npix} - |G|$ the number of its conditioning pixels. Also divide chi2 by the factor, and multiply the deviation $F z_M$ of each draw from the mean by its square root. Then, for held-out Gaussian spectra, take any fixed linear combination of the replaced pixels. Its squared residual, divided by its variance from the covariance, has the mean 1. It needs a dense, mean-centered, unweighted scaled model.
covariance (bool, optional) – If True, keep the conditional covariance of each pattern for Repaired.covariance.
rng (numpy.random.Generator, optional) – If given, the source of Repaired.draws. infill draws rng.standard_normal((n, npix)) for each chunk of $n$ rows, in row order, and nothing else. Each spectrum uses its row of these normals. Calls on consecutive slices with one generator give the draws of one call, to rounding. Draws widen the errors of measurements made on them; var and covariance give the uncertainty of the replaced pixels. None: no draws.
progress (callable, optional) – progress(stage, done, total), called after each chunk of spectra, with the stage "infill", the number of spectra done and their total. For a noise-aware model, a chunk is one spectrum.

Returns:

Repaired – The infilled spectra, in flux units.

Raises:

TypeError – If model is not a Model, or inflate or covariance is not a bool. If rng is neither None nor a numpy.random.Generator, or progress is not callable or None. If flux, ivar or mask has a dtype that scale does not allow, or replace or ignore is not of a bool dtype.
ValueError – If flux, ivar or mask has a shape or a value that scale does not allow. If flux has another npix than model, or a spectrum has no good pixel. If replace or ignore has another shape than (nspec, npix) or (npix,). If inflate is True and model has a constant center, weights, a LowRank covariance or noise. If inflate is True and $N_s - N_k - 3 \le 0$ for a spectrum with a pixel to replace, where $f_s(N_s, N_k, 1)$ is not defined; the message names $N_s$ and $N_k$ .
numpy.linalg.LinAlgError – If the covariance of a bin fails the check of its precision matrix, as in pixelwise, where a spectrum of the bin has a pixel to replace and $G$ does not cover all its pixels. If $P_{GG}$ of a pattern is not positive definite. For a pattern whose $G$ covers all pixels, if cov[b] at $T$ fails that check, or if rng is given and cov[b] at $M$ is not positive definite. The message names the bin, and for a pattern its first spectrum.

For a noise-aware model, $C_s$ replaces cov[b]. A $C_s$ that fails the check of its precision matrix raises with its row in the message. A LowRank cov[0] forms no precision matrix. Its $C_s$ raises where $\kappa_K > 10^{13}$ , with
$\kappa_K = 1 + \lambda_{\max}(\Lambda^{1/2} V^\mathsf{T} D_s^{-1} V \Lambda^{1/2}),$
$\Lambda = \mathrm{diag}(\mathrm{values})$ , and $D_s = \mathrm{diag}(\mathrm{diag} + \nu_s)$ of $C_s$ . It also raises where the factor of $I + \Lambda^{1/2} V_O^\mathsf{T} D_{s,O}^{-1} V_O \Lambda^{1/2}$ fails, or the factor of its conditional covariance at $M$ (at $T$ without rng). The message names its row. $\kappa_K$ accepts every $C_s$ that the $\kappa_W$ check of LowRank accepts. It also accepts some with $\kappa_W$ up to $\mathrm{npix} \, \kappa_K$ , which pixelwise and repair still reject.

Notes

Cost:

$\mathrm{npix}^3$ flops for $P$ at the first use of a bin, or $\mathrm{npix}^2 K$ for a LowRank of rank $K$ .
Per pattern and chunk, $|G|^3/3 + |M|^3/3$ flops, and $|M|^3/3$ more with covariance=True.
Per spectrum, $2 |G| (\mathrm{npix} + |G|)$ flops, or $2 |G| \, \mathrm{npix}$ where $G$ is one run.
Without rng, some patterns share one factor: $n \ge 2$ patterns of a chunk with one bin and one $M$ , and $G$ short of all pixels. Let $B$ be the pixels in the $G$ of all of them. For each pattern, let $E$ be the pixels of its $G$ outside $B$ . Let $D$ be the pixels of $M$ in the $T$ of another of them but not in its own. The patterns share the factor under three conditions. Each has $|E| + |D| \le |B|/16$ , and $|B| \ge 64 (1 + 2/n)$ . The shared factor costs fewer flops than their own factors, or than the mean from $C_{OO}$ .
These patterns together cost $|B|^3/3 + |M|^3/3$ flops, and $|M|^3/3$ more with covariance=True. Let $U$ be the pixels in the $T$ of any of them. Each of these patterns costs $|B|^2 |E| + |M|^2 |E| + |U| (|D| + |E|)^2$ flops, and $|M|^2 |E|$ more with covariance=True. Each of their spectra costs $2|G|\,\mathrm{npix} + |B|^2 + |M|^2 + 2|U|^2 + 4(|U| + |E|)(|D| + |E|)$ flops. An (npix,) mask on spectra with a few bad pixels of their own is such a case.
The model keeps $P$ (see Model).
With covariance=True, the result keeps $|M|^2 \times 8$ bytes per pattern and chunk.
With rng, each spectrum costs npix normals and $|M|^2$ flops more. Where $G$ covers all pixels, it costs $2 |M|^2$ flops more instead. Each pattern and chunk then also costs $|M|^3/3$ flops for the factor of cov[b] at $M$ , which a LowRank forms at $|M|^2 K$ flops.
A noise-aware model forms $C_s$ for each spectrum with a pixel to replace, and every spectrum is its own pattern. Where $G$ is also short of all pixels, a dense cov[0] forms no $P_s$ . Per spectrum, it costs $\mathrm{npix}^3/3$ flops for the Cholesky factor of $C_s$ , before the pattern costs above. From this factor, $P_{GG}$ and $P_{G,:} x_0$ cost about $(\mathrm{npix} - g_0)^2 |G| + \mathrm{npix} |G|^2$ flops, with $g_0$ the first pixel of $G$ .
A LowRank cov[0] of rank $K$ forms no $P_s$ . Per spectrum, it costs $\mathrm{npix} K^2$ flops for its $K \times K$ system, and $|M| K^2$ for the mean and var. chi2 costs $|T|^2 K + |T|^3/3$ . With rng, chi2 and the draws cost $|M|^2 K + |M|^3/3$ instead.

Repair

Detect artifacts, and replace the flagged pixels and the bad pixels.

gspice.repair

gspice.repair(
    model: Model,
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
    *,
    nguard: int | np.integer[Any] = 20,
    nsigma_replace: float | np.floating[Any] | np.integer[Any] = 10.0,
    nsigma_ignore: float | np.floating[Any] | np.integer[Any] = 5.0,
    dilate: int | np.integer[Any] = 1,
    blockwise: bool | np.bool_ = True,
    runlength: int | np.integer[Any] | None = 5,
    halfwidth: int | np.integer[Any] = 150,
    minfrac: float | np.floating[Any] | np.integer[Any] = 0.5,
    count_on: Literal['all', 'flagged'] = 'all',
    inflate: bool | np.bool_ = False,
    covariance: bool | np.bool_ = False,
    rng: np.random.Generator | None = None,
    progress: Callable[[str, int, int], object] | None = None,
) -> Repaired

Detect artifacts, and replace the flagged pixels and the bad pixels.

Take a spectrum with the working vector $x$ (see Model) and the bad pixels (mask != 0) | (ivar == 0). For a vector $y$ , let $z(y \mid w) = (y - \mu(w)) / \sigma$ at every pixel, its interpolated bad pixels included. Here $\mu(w)$ and $\sigma^2$ are the mean and variance of pixelwise with nguard for a vector $w$ , and $z(y) = z(y \mid y)$ is the z of pixelwise. Let $\mathrm{dil}(F, d)$ be the pixels within $d$ of a pixel of $F$ along the spectrum, with $d$ = dilate. Let $\tau_r$ be nsigma_replace and $\tau_i$ be nsigma_ignore. A conditional mean is that of infill. The steps, in this order, are:

$M_1 = \mathrm{dil}(|z(x)| > \tau_r, d) \cup \mathrm{bad}$ , and $I_1 = \mathrm{dil}(|z(x)| > \tau_i, d) \setminus M_1$ .
If blockwise: let $x_{\mathrm{blk}}$ be $x$ with its values at $M_1 \cup I_1$ replaced by their conditional mean given the other pixels. Take the blocks $A$ of model.blocks that are not flagged, and the ratio $\chi^2_A(x_{\mathrm{blk}}) / \mathrm{thresholds}[b, A]$ . Here $\chi^2_A$ is that of fit_blocks, and $b$ is the S/N bin of the spectrum. While the largest ratio exceeds 1, the block of that ratio is flagged, and its pixels join $M_{\mathrm{blk}}$ . $x_{\mathrm{blk}}$ then becomes $x$ with its values at $M_1 \cup I_1 \cup M_{\mathrm{blk}}$ replaced by their conditional mean given the other pixels. A block is flagged at most once, so there are at most nblock rounds. Without blockwise, $M_{\mathrm{blk}}$ is empty.
Let $x_1$ be $x$ with its values at $M_1 \cup M_{\mathrm{blk}}$ replaced by their conditional mean given the pixels outside $M_1 \cup M_{\mathrm{blk}} \cup I_1$ . Let $z_2 = z(x \mid x_1)$ . Then
$M = \mathrm{dil}\bigl((|z_2| > \tau_r) \cup ((|z(x)| > \tau_r) \cap (|z_2| > \tau_i)), d\bigr) \cup M_{\mathrm{blk}} \cup \mathrm{bad},$
and $I = \mathrm{dil}(|z_2| > \tau_i, d) \setminus M$ .
If runlength is not None, step 4 widens $I$ . Every run of at least runlength consecutive pixels of the counting set $R$ adds the pixels within halfwidth of the run to $I$ , minus $M$ . Step 4 adds no pixel to a spectrum if fewer than $\lfloor \mathrm{minfrac} \cdot \mathrm{npix} \rfloor$ pixels outside $M$ and $I$ would then remain. $R$ is $M \cup I$ for count_on="all". For "flagged", $R$ holds $M_{\mathrm{blk}}$ and the pixels outside bad with $|z_2| > \tau_i$ .
The result is that of infill with replace=M, ignore=I, and inflate, covariance and rng, bit for bit, plus the flagged blocks. A dense noise-aware model can form $P_s$ in steps 1 to 4, for guard windows or blocks short of all pixels. Step 5 then conditions on that $P_s$ , where infill uses the Cholesky factor of $C_s$ . The two then agree to rounding.

Nothing limits the replaced set. Since $M$ and $I$ depend on the data, $f_s$ does not hold for Repaired.chi2 (see finite_sample_factor). For a noise-aware model, $x$ is in flux units, and the $C_s$ of the spectrum (see Model) replaces cov[b] in every step. repair then runs the five steps one spectrum at a time, with one $P_s$ for all of them.

Parameters:

model (Model) – The model; with blockwise=True, one with blocks from fit_blocks.
flux, ivar, mask (array_like, shape (nspec, npix)) – As in scale, with the npix of model.
nguard (int, optional) – The half-width of the guard window of z (pixels); 0 or more.
nsigma_replace (float, optional) – The threshold $\tau_r$ of $|z|$ for replacing a pixel; finite and positive.
nsigma_ignore (float, optional) – The threshold $\tau_i$ of $|z|$ for ignoring a pixel; finite, positive and at most nsigma_replace.
dilate (int, optional) – The number of pixels $d$ by which dil widens the flags of z on each side; 0 or more. Paper Fig. 8 uses 2.
blockwise (bool, optional) – If True, run the block pass of step 2.
runlength (int or None, optional) – The shortest run of $R$ that widens $I$ (pixels); 1 or more. None skips step 4.
halfwidth (int, optional) – The number of pixels on each side of a run that join $I$ ; 0 or more.
minfrac (float, optional) – The fraction of npix that step 4 leaves to condition on; from 0 to 1.
count_on ({“all”, “flagged”}, optional) – The counting set $R$ of step 4.
inflate, covariance (bool, optional) – As in infill, for step 5. inflate needs a dense, mean-centered, unweighted scaled model.
rng (numpy.random.Generator, optional) – As in infill, for step 5: the source of Repaired.draws. None: no draws.
progress (callable, optional) – progress(stage, done, total), called with the stage "repair", the number of spectra done and their total. repair calls it after each chunk of spectra of steps 1 to 4, and again after each chunk of step 5. A noise-aware model runs the five steps on chunks of one spectrum, and reports after its steps 1 to 4 and after its step 5.

Returns:

Repaired – The repaired spectra, in flux units, with replaced $M$ and ignored $I$ . blocks holds the flagged blocks of each spectrum, or None if blockwise is False.

Raises:

TypeError – If model is not a Model, or nguard, dilate, runlength or halfwidth is not an integer (or None for runlength). If nsigma_replace, nsigma_ignore or minfrac is not a real number, or blockwise, inflate or covariance is not a bool. If rng is neither None nor a numpy.random.Generator, or progress is not callable or None. If flux, ivar or mask has a dtype that scale does not allow.
ValueError – If flux, ivar or mask has a shape or a value that scale does not allow. If flux has another npix than model, or a spectrum has no good pixel. If an argument is outside its range above, count_on is neither "all" nor "flagged", or blockwise is True and model has no blocks. If inflate is True and model has a constant center, weights, a LowRank covariance or noise. If inflate is True and, before step 5, $N_s - N_k - 3 \le 0$ for a spectrum with a replaced pixel, with $N_k = \mathrm{npix} - |M \cup I|$ ; the message names $N_s$ and $N_k$ .
numpy.linalg.LinAlgError – If the covariance of a bin with spectra fails the check of its precision matrix, as in pixelwise. repair checks the $P$ of every bin with spectra before any spectrum, also for nguard >= npix - 1, where pixelwise computes no $P$ . As in fit_blocks for the blocks of a bin with spectra, with blockwise=True. As in infill for a conditional mean. The message names the bin, and the block or the row. For a noise-aware model, as in infill, with $C_s$ in place of cov[b], checked for each spectrum where a step factors $C_s$ .

Notes

Cost:

$\mathrm{npix}^3$ flops for $P$ at the first use of a bin, and the guard weights as in pixelwise.
With blockwise=True, $|G|^3/3$ flops per block and bin with spectra, with $G$ the block and its guard.
Per spectrum, $4 \, \mathrm{npix}^2$ flops for the two pixelwise passes, and about $2 \, \mathrm{npix}^2 + \sum_A |G|^2$ flops per scan of its blocks.
Per spectrum, one scan and one infill per flagged block, and two infills as in infill. These two are the infill of step 5, and that of step 2, whose last $x_{\mathrm{blk}}$ gives $x_1$ , or without blockwise that of step 3. All infills except that of step 5 compute only the mean, from $C_{OO}$ where $|G| > \mathrm{npix}/2$ .
repair holds replaced and ignored of all spectra, and a few arrays of one chunk.
A noise-aware model costs, for each spectrum, $P_s$ , its guard windows as in pixelwise, and its block factors. $P_s$ takes $\mathrm{npix}^3$ flops for a dense cov[0], or $\mathrm{npix}^2 K$ for a LowRank of rank $K$ . The infills of a LowRank cov[0] do not use $P_s$ , and cost as in infill.

Diagnostics

Check that the predicted uncertainties are calibrated.

gspice.Calibration

class gspice.Calibration(z2: ArrayLike, n: ArrayLike, fs: ArrayLike)

Held-out whitened residuals, from calibration or another held-out result.

Row $s$ is a spectrum, and column $j$ is one prediction of n[s, j] of its pixels. z2[s, j] is $r^\mathsf{T} C^{-1} r$ , with $r$ the residuals of these pixels from their conditional mean and $C$ their full conditional covariance. fs[s, j] is its finite-sample factor, such as $f_s$ of finite_sample_factor. For held-out Gaussian spectra and pixel sets that do not depend on the data, z2 / (n * fs) then has the mean 1. The methods use only the entries with n > 0 and finite z2 and fs.

Calibration(z2=r.chi2[:, None], n=r.dof[:, None], fs=r.fs[:, None]) holds the result r of infill on held-out spectra, with replace and ignore masks that do not depend on the data. Calibration(z2=p.z**2, n=1, fs=1) holds the result p of pixelwise with inflate=True.

The constructor takes array_like values, and broadcasts n and fs to the shape of z2. It stores read-only copies: z2 and fs as float64, and n as int64. Two calibrations are equal only if they are the same object. Pickle, copy and deepcopy check the fields again.

Attributes:

z2 (ndarray, shape (nspec, ncol)) – $Z^2$ of each prediction (unit-free): finite and 0 or more, or NaN.
n (ndarray of int64, shape (nspec, ncol)) – The number of predicted pixels: from 0 to $2^{53} - 1$ .
fs (ndarray, shape (nspec, ncol)) – The finite-sample factor of each prediction (unit-free): finite and positive, or NaN.

Raises:

TypeError – If z2 or fs does not have an integer or floating dtype, or n does not have an integer dtype; bool is neither. The message starts with the name of the field.
ValueError – If z2 is not 2-D, or n or fs does not broadcast to the shape of z2. If a field has a value that is not allowed above. The message starts with the name of the field.

ratio

ratio(where: ArrayLike | None = None) -> tuple[NDArray[np.float64], NDArray[np.float64]]

Mean of z2 / (n * fs) in each column, and its standard error.

Take the $m$ entries of column $j$ that the methods use, in the spectra of where. mean[j] is the mean of their z2 / (n * fs), and $\mathrm{sem}[j] = \mathrm{std} / \sqrt{m}$ , with std their standard deviation with the divisor $m - 1$ . For held-out spectra, sem is conditional on the fitted model. It falls with the number of held-out spectra, but the spread across training sets does not.

Parameters:

where (array_like of bool, shape (nspec,), optional) – The spectra to use. None: all.

Returns:

mean (ndarray, shape (ncol,)) – The mean (unit-free); NaN where $m = 0$ .
sem (ndarray, shape (ncol,)) – Its standard error (unit-free); NaN where $m < 2$ .

Raises:

TypeError – If where does not have a bool dtype.
ValueError – If where has another shape than (nspec,).

tail

tail(
    q: float | np.floating[Any] | np.integer[Any] = 0.99,
    where: ArrayLike | None = None,
) -> NDArray[np.float64]

Fraction of each column with z2 / fs above a $\chi^2_n$ quantile.

Take the $m$ entries of column $j$ that the methods use, in the spectra of where. The result is the fraction of them with z2 / fs above $F_n^{-1}(q)$ . $F_n$ is the $\chi^2$ distribution function with the n of the entry as its degrees of freedom. For held-out Gaussian spectra, the fraction approaches $1 - q$ only for $N_s - N_k \gg n$ . It is 1.37% at $(N_s, N_k, n) = (200, 40, 10)$ for $q = 0.99$ .

Parameters:

q (float, optional) – The quantile; $0 < q < 1$ .
where (array_like of bool, shape (nspec,), optional) – The spectra to use. None: all.

Returns:

ndarray, shape (ncol,) – The fraction (unit-free); NaN where $m = 0$ .

Raises:

TypeError – If q is not a real number, or where does not have a bool dtype.
ValueError – If q is not strictly between 0 and 1, or where has another shape than (nspec,).

qq

qq(
    i: int | np.integer[Any] = 0,
    where: ArrayLike | None = None,
) -> tuple[NDArray[np.float64], NDArray[np.float64]]

Normal scores of column i against $\chi^2_n$ , for a QQ plot.

Take the $m$ entries of column i that the methods use, in the spectra of where. expected holds $\Phi^{-1}((j - 1/2) / m)$ for $j = 1, \dots, m$ , with $\Phi$ the standard normal distribution function. observed holds the sorted scores $\Phi^{-1}(F_n(\mathrm{z2} / \mathrm{fs}))$ . $F_n$ is the $\chi^2$ distribution function with the n of the entry as its degrees of freedom. A score is $\Phi^{-1}(F_n)$ where $F_n < 1/2$ . It is $-\Phi^{-1}(1 - F_n)$ elsewhere, from the survival function, so that both tails keep their precision. A score is $-\infty$ where z2 is 0, and $+\infty$ where the survival function underflows to 0. For held-out Gaussian spectra, observed approaches expected under the condition of tail.

Parameters:

i (int, optional) – The column; from 0 to ncol − 1.
where (array_like of bool, shape (nspec,), optional) – The spectra to use. None: all.

Returns:

expected (ndarray, shape (m,)) – The normal scores of the plotting positions, in increasing order.
observed (ndarray, shape (m,)) – The scores of the entries, in increasing order.

Raises:

TypeError – If i is not an integer, or where does not have a bool dtype.
ValueError – If i is outside 0 to ncol − 1, or where has another shape than (nspec,).

gspice.calibration

gspice.calibration(
    model: Model,
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
    *,
    sizes: Iterable[int | np.integer[Any]] | None = None,
    layout: Literal['run', 'scattered'] = 'run',
    rng: np.random.Generator,
    progress: Callable[[str, int, int], object] | None = None,
) -> Calibration

Held-out whitened residuals for masks drawn independently of the data.

For each spectrum and each size $n$ of sizes, calibration draws a set $M$ of $n$ good pixels. It computes, as infill with replace=M does,

z^2 = r_M^\mathsf{T} C_{MM}^{-1} r_M,

with $r_M$ the input values at $M$ minus their conditional mean given the pixels outside $G = M \cup \mathrm{bad}$ , and $C_{MM}$ their full conditional covariance. It also computes $f_s(N_s, \mathrm{npix} - |G|, n)$ of finite_sample_factor, NaN where it is not defined. $N_s$ is the n_train of the bin of the spectrum. The factor is 1 where no factor is derived (see Repaired.fs). For held-out Gaussian spectra of a model with a derived factor, $z^2 / (n f_s)$ then has the mean 1. For a noise-aware model, the covariance is the $C_s$ of the spectrum (see Model).

Rank the $m$ good pixels of a spectrum from 0 to $m - 1$ in increasing order. With layout="run", $M$ holds those of ranks start to start + n - 1, with start = rng.integers(m - n + 1). With layout="scattered", $M$ holds those of ranks rng.choice(m, n, replace=False). rng draws these spectrum by spectrum in row order, and size by size within a spectrum. It draws nothing else, so the masks do not depend on the chunking. A spectrum with $m < n$ draws nothing for that size. Its entry there has n 0 and z2 NaN.

Parameters:

model (Model) – The model.
flux, ivar, mask (array_like, shape (nspec, npix)) – As in scale, with the npix of model.
sizes (iterable of int, optional) – The sizes $n$ , one column of the result each, in this order; 1 or more. None: the sizes 1, 4, 16, 64, 256 and 1024 that are not above half the median number of good pixels per spectrum (np.median). Without a spectrum, all six. A spectrum with fewer good pixels than a size has n 0 there, so the spectra of a column can differ between sizes. Calibration methods with where=(result.n > 0).all(axis=1) compare all sizes on the same spectra.
layout ({“run”, “scattered”}, optional) – Whether $M$ is a run of good pixels, or good pixels scattered at random.
rng (numpy.random.Generator) – The source of the masks, as above.
progress (callable, optional) – progress(stage, done, total), called after each chunk of spectra, with the stage "calibration", the number of spectra done and their total. For a noise-aware model, a chunk is one spectrum.

Returns:

Calibration – z2, n and fs of each spectrum (row) and size (column).

Raises:

TypeError – If model is not a Model, or sizes is neither None nor an iterable of integers. If rng is not a numpy.random.Generator, or progress is not callable or None. If flux, ivar or mask has a dtype that scale does not allow.
ValueError – If flux, ivar or mask has a shape or a value that scale does not allow. If flux has another npix than model, or a spectrum has no good pixel. If sizes is empty or holds an integer below 1, or layout is neither "run" nor "scattered". If sizes is None and the median number of good pixels per spectrum is below 2.
numpy.linalg.LinAlgError – Before the first chunk, if the covariance of a bin fails its check, for a bin with a spectrum that has more good pixels than the smallest size. For a noise-aware model, the check is at the spectrum. Otherwise as in infill with replace=M.

Notes

Cost:

For each size, $\mathrm{npix}^3$ flops for $P$ at the first use of a bin.
Then $|G|^3/3$ flops per pattern of a chunk of max(1, 2**24 // npix) spectra, or $n^3/3$ where $G$ covers all pixels. A spectrum has at most one pattern per size.
Per spectrum, $2 |G| (\mathrm{npix} + |G|)$ flops, or $2 |G| \, \mathrm{npix}$ where $G$ is one run.
The masks of a chunk take up to $\mathrm{len(sizes)} \times 2^{24}$ bytes.
A noise-aware model costs per spectrum and size as in infill, with $M = T$ . A dense cov[0] factors $C_s$ once per spectrum for all sizes.

gspice.finite_sample_factor

gspice.finite_sample_factor(
    n_train: int | np.integer[Any],
    n_cond: int | np.integer[Any],
    n_pred: int | np.integer[Any] = 1,
) -> float
gspice.finite_sample_factor(
    n_train: ArrayLike,
    n_cond: ArrayLike,
    n_pred: ArrayLike = 1,
) -> float | NDArray[np.float64]

Held-out calibration factor $E[Z^2] / N_*$ for Gaussian data.

Take a held-out spectrum, predicted at $N_*$ pixels from $N_k$ conditioning pixels. The prediction uses the mean and the covariance (divisor $N_s - 1$ ) of $N_s$ unweighted training spectra. Let $r$ be the residual of the $N_*$ pixels from their conditional mean, and $C$ their full conditional covariance. Then $Z^2 = r^\mathsf{T} C^{-1} r$ has the mean

\frac{E[Z^2]}{N_*} = f_s = \left(1 + \frac{1}{N_s}\right) \frac{N_s - 2}{N_s - N_k - 2} \, \frac{N_s - 1}{N_s - N_k - N_* - 2}.

The factor is exact for independent Gaussian working vectors and pixel sets that do not depend on the data. Pixels that are neither conditioned on nor predicted, such as a guard window or ignored pixels, do not count in $N_k$ . $f_s$ does not hold for data that are not Gaussian, such as scaled spectra with a spread of S/N. The arguments broadcast. The arithmetic is float64, which holds each count exactly.

Parameters:

n_train (array_like of int) – $N_s$ , the number of training spectra; from 1 to $2^{53} - 1$ .
n_cond (array_like of int) – $N_k$ , the number of conditioning pixels; from 0 to $2^{53} - 1$ .
n_pred (array_like of int, optional) – $N_*$ , the number of predicted pixels; from 1 to $2^{53} - 1$ .

Returns:

float or ndarray – $f_s$ (unit-free): a float if the arguments are scalars, and otherwise a float64 array of their broadcast shape.

Raises:

TypeError – If an argument does not have an integer dtype; bool is not one.
ValueError – If an argument is outside its range, the arguments do not broadcast, or $N_s - N_k - N_* - 2 \le 0$ . The message names the arguments and the first such entry.

gspice.noise_audit

gspice.noise_audit(
    model: Model,
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
    *,
    nguard: int | np.integer[Any] = 20,
    progress: Callable[[str, int, int], object] | None = None,
) -> NDArray[np.float64]

Ratio of the residual scatter to the reported noise at each pixel.

Let $f$ be the flux and $\hat f$ its prediction by pixelwise with nguard, in flux units. Then

\rho_i^2 = \operatorname{mean}_s \, (f_{si} - \hat f_{si})^2 \, \mathrm{ivar}_{si},

where the mean runs over the spectra $s$ whose pixel $i$ is good (see scale). Let $t$ be the noise-free flux. Let $\eta_i$ be the mean of $\mathrm{ivar}_{si} \operatorname{var}(f_{si} - t_{si})$ , which is 1 where ivar is right, and $e_i$ the mean of $\mathrm{ivar}_{si} (t_{si} - \hat f_{si})^2$ . Then

E[\rho_i^2] = \eta_i + E[e_i]

if $\hat f_{si}$ does not depend on the noise of $f_{si}$ . That holds for held-out spectra, with mask and ivar independent of the noise, and noise that is uncorrelated beyond nguard. It fails where the interpolation of a bad pixel carries $f_{si}$ , or a pixel within nguard of $i$ , into $\hat f_{si}$ (see pixelwise). With several S/N bins, $f_{si}$ also enters the median S/N that picks the bin of $s$ .

$e_i$ holds the astrophysical error, the noise propagated from the conditioning pixels, and the estimation error of model. So $E[\rho_i^2] \ge \eta_i$ where the identity holds, and $e_i$ grows with S/N. Next to a run of bad pixels longer than nguard, $\rho_i^2$ can fall below $\eta_i$ . A scaled model absorbs into its covariance an $\eta_i$ that its training spectra share, so the z of pixelwise is blind to it. $\rho$ is not.

Parameters:

model (Model) – The model.
flux, ivar, mask (array_like, shape (nspec, npix)) – As in scale, with the npix of model.
nguard (int, optional) – The half-width of the guard window of pixelwise (pixels); 0 or more.
progress (callable, optional) – progress(stage, done, total), called after each chunk of spectra, with the stage "noise_audit", the number of spectra done and their total. For a noise-aware model, a chunk is one spectrum.

Returns:

ndarray, shape (npix,) – $\rho$ (unit-free); NaN at the pixels where no spectrum is good.

Raises:

TypeError – If model is not a Model, nguard is not an integer, or progress is not callable or None. If flux, ivar or mask has a dtype that scale does not allow.
ValueError – If flux, ivar or mask has a shape or a value that scale does not allow. If flux has another npix than model, nguard is negative, or a spectrum has no good pixel.
numpy.linalg.LinAlgError – As in pixelwise.

Notes

Cost: that of pixelwise. Beyond the inputs, and their float64 copies if they are long double, noise_audit holds a bool array of the input shape. It also holds a few float64 arrays of one chunk of max(1, 2**24 // npix) spectra.

Scaling

Scale spectra before fitting, and undo the scaling.

gspice.scale

gspice.scale(
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
    *,
    center: float | np.floating[Any] | np.integer[Any] = 0.0,
    refmean: ArrayLike | None = None,
) -> tuple[NDArray[np.float64], NDArray[np.float64], NDArray[np.float64]]

Scale spectra to unit mean noise.

X_s = a_s (\tilde f_s - \mathrm{center}) - \mathrm{refmean}

The bad pixels are (mask != 0) | (ivar == 0), and scale never reads flux or ivar there. It checks and computes the float64 values of the inputs. A long double beyond the float64 range is not finite, and one that rounds to 0 in float64 is 0. $\tilde f_s$ is flux[s] cast to float64, then linearly interpolated over its bad pixels, and held constant beyond its first and last good pixels. $a_s$ is the square root of the mean of $\mathrm{ivar}_s$ over its good pixels. With float64 flux, scale(flux, ivar, mask) equals the 2022 gspice.legacy.standard_scale with the mask (mask != 0) | (ivar == 0), to rounding.

Parameters:

flux (array_like, shape (nspec, npix)) – Flux, of an integer or floating dtype; finite at good pixels.
ivar (array_like, shape (nspec, npix)) – Inverse variance (flux⁻²), of an integer or floating dtype; finite and 0 or more at pixels not in mask.
mask (array_like, shape (nspec, npix), optional) – A nonzero value marks a pixel in mask; of a bool, integer or floating dtype. None: no pixel is in mask.
center (float, optional) – A constant subtracted before scaling (flux).
refmean (array_like, shape (npix,), optional) – The mean subtracted after scaling (scaled); finite. None: the column mean of $a_s (\tilde f_s - \mathrm{center})$ , which needs at least one spectrum. For held-out spectra, pass the model mean of their S/N bin: their own mean breaks the finite-sample factor $f_s$ .

Returns:

X (ndarray, shape (nspec, npix)) – The scaled spectra (scaled).
refscale (ndarray, shape (nspec,)) – $a_s$ (flux⁻¹).
refmean (ndarray, shape (npix,)) – The subtracted mean (scaled); a copy if refmean is given.

Raises:

TypeError – If an argument has a dtype or a type that is not allowed above.
ValueError – If an argument has a shape or a value that is not allowed above, or a spectrum has no good pixel.

gspice.unscale

gspice.unscale(
    x: ArrayLike,
    refscale: ArrayLike,
    refmean: ArrayLike,
    *,
    center: float | np.floating[Any] | np.integer[Any] = 0.0,
    var: ArrayLike | None = None,
) -> tuple[NDArray[np.float64], NDArray[np.float64] | None]

Invert scale, for values and their variances.

f = \frac{x + \mathrm{refmean}}{\mathrm{refscale}} + \mathrm{center}, \qquad \mathrm{var}_f = \frac{\mathrm{var}}{\mathrm{refscale}^2}

It computes with the float64 values of the inputs, as scale does.

Parameters:

x (array_like, shape (nspec, n)) – Values at $n$ pixels (scaled), of an integer or floating dtype.
refscale (array_like, shape (nspec,)) – $a_s$ of scale (flux⁻¹); finite and positive.
refmean (array_like, shape (n,)) – The mean that scale subtracted at these pixels (scaled); finite.
center (float, optional) – The center of scale (flux).
var (array_like, shape (n,) or (nspec, n), optional) – Variances of x (scaled²), of an integer or floating dtype; 0 or more.

Returns:

f (ndarray, shape (nspec, n)) – Flux.
var_f (ndarray, shape (nspec, n), or None) – Variances of f (flux²); None if var is None.

Raises:

TypeError – If an argument has a dtype or a type that is not allowed above.
ValueError – If an argument has a shape or a value that is not allowed above.

gspice.snr

gspice.snr(
    flux: ArrayLike,
    ivar: ArrayLike,
    mask: ArrayLike | None = None,
) -> NDArray[np.float64]

Median signal-to-noise ratio over the good pixels of each spectrum.

The ratio at a pixel is $\mathrm{flux} \sqrt{\mathrm{ivar}}$ . The bad pixels are (mask != 0) | (ivar == 0), and snr never reads flux or ivar there. The products are computed in float64.

Parameters:

flux, ivar, mask (array_like, shape (nspec, npix)) – As in scale.

Returns:

ndarray, shape (nspec,) – Signal-to-noise ratio (unit-free).

Raises:

TypeError – If an argument has a dtype that scale does not allow.
ValueError – If an argument has a shape or a value that scale does not allow, or a spectrum has no good pixel.

Simulation

Make synthetic spectra for tests and examples.

gspice.Simulation

class gspice.Simulation(
    flux: NDArray[np.float64],
    ivar: NDArray[np.float64],
    truth: NDArray[np.float64],
    artifacts: NDArray[np.bool_],
    snr: NDArray[np.float64],
    labels: NDArray[np.int64],
)

Synthetic spectra from simulate, with their noise-free flux and artifacts.

The arrays that simulate returns are read-only.

Attributes:

flux (ndarray, shape (nspec, npix)) – The flux, with noise and artifacts (flux).
ivar (ndarray, shape (nspec, npix)) – The inverse variance of the noise, $1/\sigma^2$ (flux⁻²).
truth (ndarray, shape (nspec, npix)) – The noise-free flux without artifacts (flux).
artifacts (ndarray of bool, shape (nspec, npix)) – The pixels that an injected artifact covers.
snr (ndarray, shape (nspec,)) – The continuum signal-to-noise ratio that each spectrum was drawn with (unit-free).
labels (ndarray of int64, shape (nspec,)) – The population of each spectrum: 0, or 1 for the giants.

gspice.simulate

gspice.simulate(
    nspec: int | np.integer[Any],
    npix: int | np.integer[Any] = 500,
    *,
    rng: np.random.Generator,
    artifacts: float | np.floating[Any] | np.integer[Any] = 0.0,
    gaussian: bool | np.bool_ = False,
) -> Simulation

Synthetic continuum-normalized spectra for tests and examples (no network).

Each spectrum has the continuum 1 and $\mathrm{nline} = \lceil \mathrm{npix}/25 \rceil$ absorption lines $j = 1, \dots, \mathrm{nline}$ at evenly spaced centers. Line $j$ has a Gaussian profile $p_j$ , with a standard deviation from 1.5 to 4 pixels. The lines depend only on npix. Two latent parameters $\theta$ set the depths. A spectrum is in population 1, the giants, with the probability 0.3, and then has $\theta \sim N((0, 2), I)$ . Otherwise $\theta \sim N(0, I)$ . Line $j$ has the depth

\mathrm{depth}_j = 0.95\,\mathrm{expit}\bigl(\mathrm{logit}(d_j / 0.95) + r_j \cdot \theta\bigr),

with $d_j$ from 0.1 to 0.6 its depth at $\theta = 0$ , and $r_j = 0.5 (\cos j, \sin j)$ . The depths are thus correlated, in (0, 0.95), and not Gaussian. Then $\mathrm{truth} = 1 - \sum_j \mathrm{depth}_j\, p_j$ , in (0, 1]. The continuum S/N has the density $p(\mathrm{S/N}) \propto (\mathrm{S/N})^{-2.5}$ on [30, 400]. The noise is Gaussian with $\sigma = \sqrt{\mathrm{truth}} / \mathrm{snr}$ .

With gaussian=True, $\mathrm{truth} = \mu + J \theta$ with $\theta \sim N(0, I)$ . Here $\mu$ is the truth of population 0 at $\theta = 0$ , and the columns of $J$ are its derivatives there. Every spectrum has the label 0 and the S/N 50, and $\sigma = \sqrt{\mu} / 50$ . So ivar is equal for all spectra. Without artifacts, $\mathrm{flux} \sim N(\mu, J J^\mathsf{T} + \mathrm{diag}(\sigma^2))$ , and the scaled spectra (see scale) are exactly Gaussian. truth can then leave (0, 1].

round(artifacts * nspec) spectra (half to even), chosen at random, get one artifact each. Its type is drawn uniformly from three:

A cosmic ray adds $h \sigma$ to each of 1 to 3 adjacent pixels, with $h$ log-uniform on [5, 50].
A bad column sets the flux of 5 to 20 adjacent pixels to 2.
A shifted segment of 200 to 600 adjacent pixels $i$ has the flux $\mathrm{truth}(i - d)$ plus the noise. The shift $d$ is from 0.5 to 3 pixels, of either sign. truth is interpolated linearly, and held constant beyond its ends.

An artifact covers at most npix − 1 pixels, or 1 if npix is 1. Its start is drawn uniformly from the starts where it fits. For one state of rng, every value of artifacts gives the same ivar, truth, snr and labels. It also gives the same flux outside the pixels in Simulation.artifacts.

Parameters:

nspec (int) – The number of spectra; 0 or more.
npix (int, optional) – The number of pixels; 1 or more.
rng (numpy.random.Generator) – The source of all randomness. With the same libraries and thread settings, the same state gives the same result bit for bit.
artifacts (float, optional) – The fraction of spectra with an artifact; from 0 to 1.
gaussian (bool, optional) – Draw exactly Gaussian spectra, as above.

Returns:

Simulation – The spectra.

Raises:

TypeError – If an argument has a type that is not allowed above.
ValueError – If an argument has a value that is not allowed above.

Notes

Cost: $O(\mathrm{nspec}\,\mathrm{npix})$ time and memory, plus the nline profiles of npix values each. With gaussian=False, the line sum adds $2\,\mathrm{nspec}\,\mathrm{nline}\,\mathrm{npix}$ flops.

Input and output

Write a model to a file, and read it back.

gspice.save

gspice.save(
    model: Model,
    path: str | bytes | os.PathLike[str] | os.PathLike[bytes],
) -> None

Write model to the .npz file path.

The archive is not compressed. It holds these arrays:

schema (the integer 1), gspice_version, mean, n_train, edges, center (NaN if None) and weighted. For a noise-aware model, also noise (True).
For the covariance of each bin $b$ , cov{b}, or vectors{b}, values{b} and diag{b} for the fields of a LowRank.
If model.blocks is not None, block_starts, block_width, block_nguard and block_thresholds for its fields.

load reads it.

Parameters:

model (Model) – The model to write.
path (str or path-like) – The file to write, with any suffix: save adds none. An existing file is replaced.

Raises:

TypeError – If model is not a Model, or path is not a str or path-like.
OSError – If the file cannot be written.

gspice.load

gspice.load(path: str | bytes | os.PathLike[str] | os.PathLike[bytes]) -> Model

Read a model written by save.

Parameters:

path (str or path-like) – A file written by save.

Returns:

Model – The model.

Raises:

ValueError – If the file is not a zip archive, or is a damaged one. If it has no schema, or a schema other than the 0-d integer 1. If it lacks one of the arrays schema, mean, n_train, edges, center and weighted. If it holds some but not all of the four block arrays, or vectors{b} without values{b} or diag{b}. If it holds an array that save does not write. The message names path. Also if Model, LowRank or Blocks rejects a field, for example a count of covariances other than the number of bins. The message then names the field.
TypeError – If path is not a str or path-like, or if Model, LowRank or Blocks rejects a field.
OSError – If the file cannot be read, and for some damaged archives.