pysersic

Submodules

Classes

`PySersicMultiPrior`	Class used for priors for multi source fitting in PySersic
`PySersicSourcePrior`	Class used for priors for single source fitting in PySersic
`SourceProperties`	A Class used to estimate initial guesses for source properties.
`FitMulti`	Class used to fit multiple sources within a single image
`FitSingle`	Class used to fit a single source
`FourierRenderer`	Class to render sources based on rendering them in Fourier space. Sersic profiles are modeled as a series of Gaussian following Shajib (2019) (https://arxiv.org/abs/1906.08263) and the implementation in lenstronomy (https://github.com/lenstronomy/lenstronomy/blob/main/lenstronomy/LensModel/Profiles/gauss_decomposition.py)
`HybridRenderer`	Class to render sources based on the hybrid rendering scheme introduced in Lang (2020). This avoids some of the artifacts introduced by rendering sources purely in Fourier space. Sersic profiles are modeled as a series of Gaussian following Shajib (2019) (https://arxiv.org/abs/1906.08263) and the implementation in lenstronomy (https://github.com/lenstronomy/lenstronomy/blob/main/lenstronomy/LensModel/Profiles/gauss_decomposition.py).
`PixelRenderer`	Render class based on rendering in pixel space and then convolving with the PSF
`PySersicResults`

Functions

`autoprior`(→ PySersicSourcePrior)	Simple wrapper function to generate a prior using the built-in defaults. This can be used as a starting place but may not work for all sources
`check_input_data`(data, rms, psf[, mask])	Check input data for certain conditions and raise warnings or exceptions if needed
`parse_mask`([mask, data])
`parse_multi_results`(→ PySersicResults)	Function written to parse results from a FitMulti instance to isolate a single source. A new PySersicResults class is created with only the posteriors of the specified source. The original chains saved under .idata_all

Package Contents

Bases: BasePrior

Class used for priors for multi source fitting in PySersic

catalog: dict

N_sources: int

__repr__() → str

class pysersic.PySersicSourcePrior(profile_type: str, sky_type: str | None = 'none', sky_guess: float | None = None, sky_guess_err: float | None = None, suffix: str | None = '')

Bases: BasePrior

Class used for priors for single source fitting in PySersic

profile_type: str

__repr__() → str

check_vars(verbose=False) → bool

Function to check if all variable for the specified profile type are set with no extras

Parameters:: verbose (bool, optional) – Wheter to print out missing and extra variables, by default False
Returns:: True if all variable for given type are present with no extra, False otherwise
Return type:: bool

class pysersic.SourceProperties(image: numpy.array | jax.numpy.array, mask: numpy.array | jax.numpy.array = None)

A Class used to estimate initial guesses for source properties. If no guesses are provided, then the class will estimate them using the photutls package and the data_properties() function.

measure_properties(**kwargs) → SourceProperties

Measure default properties of the source

Returns:: returns self
Return type:: SourceProperties

set_sky_guess(sky_guess: float | None = None, sky_guess_err: float | None = None, n_pix_sample: int = 5, **kwargs) → SourceProperties

Measure or set guess for initial sky background level. If no estimate is provided, the median of the n_pix_sample number of pixels around each edge is used

Parameters:

sky_guess (Optional[float], optional) – Initial guess for level of background, by default None
sky_guess_err (Optional[float], optional) – Uncertainity on inital guess, by default None
n_pix_sample (int, optional) – Number of pixels around each edge to use to estimate sky level if neccesary, by default 5

Returns:

returns self

Return type:

SourceProperties

set_flux_guess(flux_guess: float | None = None, flux_guess_err: float | None = None, **kwargs) → SourceProperties

Measure or set guess for initial flux. If no estimate is provided, the flux of the source in estimated as the total flux within the sgmentatated region for the source

Parameters:

flux_guess (Optional[float], optional) – Initial guess for flux, by default None
flux_guess_err (Optional[float], optional) – Uncertainty on initial guess, by default None

Returns:

returns self

Return type:

SourceProperties

set_r_eff_guess(r_eff_guess: float | None = None, r_eff_guess_err: float | None = None, **kwargs) → SourceProperties

Measure or set guess for effective radius. If no estimate is provided, the r_eff of the source in estimated using photutils

Parameters:

r_eff_guess (Optional[float], optional) – Initial guess for effective radius, by default None
r_eff_guess_err (Optional[float], optional) – Uncertainty on initial guess, by default None

Returns:

returns self

Return type:

SourceProperties

set_theta_guess(theta_guess: float | None = None, **kwargs) → SourceProperties

Measure or set guess for initial position angle. If no estimate is provided, the position angle of the source in estimated using the data_properties() function from photutils

Parameters:: theta_guess (Optional[float], optional) – Estimate of the position angle in radians, by default None
Returns:: returns self
Return type:: SourceProperties

set_position_guess(position_guess: Iterable[float, float] | None = None, **kwargs) → SourceProperties

Measure or set guess for initial position. If no estimate is provided, the position of the source in estimated using the data_properties() function from photutils

Parameters:: position_guess (Optional[Iterable[float,float]], optional) – A 2 element list, tuple or array which contain the x,y pixel values of the inital guess for the centroid, by default None
Returns:: returns self
Return type:: SourceProperties

generate_prior(profile_type: str, sky_type: str | None = 'none', suffix: str | None = '') → PySersicSourcePrior

Function to generate default priors based on a given image and profile type

Parameters:

profile_type (str) – Type of profile
sky_type (str, default 'none') – Type of sky model to use, must be one of: ‘none’, ‘flat’, ‘tilted-plane’
suffix (str, default '') – Add suffix onto all source related variables, generally only used to number sources in MultiPrior

Returns:

Pysersic prior object to be used to initialize FitSingle

Return type:

PySersicSourcePrior

visualize(figsize: Tuple[float, float] = (6.0, 6.0), cmap: str = 'gray', scale: float = 1.0) → None

Display a figure summarizing the current guess for the source properties

Parameters:

figsize (Tuple[float,float], optional) – figure, by default (6.,6.)
cmap (str, optional) – color map, by default ‘gray’
scale (float, optional) – number of +/- std’s at which to clip image, by default 1

pysersic.autoprior(image: numpy.array, profile_type: str, mask: numpy.array = None, sky_type: str = 'none') → PySersicSourcePrior

Simple wrapper function to generate a prior using the built-in defaults. This can be used as a starting place but may not work for all sources

Parameters:

image (np.array) – science image
profile_type (str) – Type of profile to be fit
mask (np.array, optional) – pixel by pixel mask, by default None
sky_type (str) – Type of sky to fit, default ‘none’

Returns:

Prior object that can be used in initializing FitSingle

Return type:

PySersicSourcePrior

class pysersic.FitMulti(data: jax.typing.ArrayLike, rms: jax.typing.ArrayLike, psf: jax.typing.ArrayLike, prior: pysersic.priors.PySersicMultiPrior, mask: jax.typing.ArrayLike | None = None, loss_func: Callable | None = gaussian_loss, renderer: pysersic.rendering.BaseRenderer | None = HybridRenderer, renderer_kwargs: dict | None = {})

Bases: BaseFitter

Class used to fit multiple sources within a single image

build_model(return_model: bool = True) → Callable

Generate Numpyro model for the specified image, profile and priors

Returns:: model – Function specifying the current model in Numpyro, can be passed to inference algorithms
Return type:: Callable

find_MAP(return_model: bool | None = True, rkey: jax.random.PRNGKey | None = jax.random.PRNGKey(3), purge_extra: bool | None = True)

Find the “best-fit” parameters as the maximum a-posteriori and return a dictionary with values for the parameters.

Parameters:

return_model (Optional[bool], optional) – whether to return the model image, adds a small time and memory overhead, by default True
rkey (Optional[jax.random.PRNGKey], optional) – rng key, by default jax.random.PRNGKey(3)

Returns:

dictionary with fit parameters and their values.

Return type:

dict

class pysersic.FitSingle(data: jax.typing.ArrayLike, rms: jax.typing.ArrayLike, psf: jax.typing.ArrayLike, prior: pysersic.priors.PySersicSourcePrior, mask: jax.typing.ArrayLike | None = None, loss_func: Callable | None = gaussian_loss, renderer: pysersic.rendering.BaseRenderer | None = HybridRenderer, renderer_kwargs: dict | None = {})

Bases: BaseFitter

Class used to fit a single source

build_model(return_model: bool = True) → Callable

Generate Numpyro model for the specified image, profile and priors

Returns:: model – Function specifying the current model in Numpyro, can be passed to inference algorithms
Return type:: Callable

pysersic.check_input_data(data: jax.typing.ArrayLike, rms: jax.typing.ArrayLike, psf: jax.typing.ArrayLike, mask: jax.typing.ArrayLike = None)

Check input data for certain conditions and raise warnings or exceptions if needed

Parameters:

data (ArrayLike) – input image (galaxy/cutout)
rms (ArrayLike) – rms/error map of the data.
psf (ArrayLike) – pixelized PSF
mask (ArrayLike, optional) – mask with True/1 indicating a pixel should be masked, by default None

Raises:

RMSWarning – If the rms map is highly discrepant from the rms of the data
PSFNormalizationWarning – if the pst normalization is not ~1
KernelError – if the provided PSF is larger than the input image (exception)
MaskWarning – If more than 70% of the image is masked.

pysersic.parse_mask(mask: jax.typing.ArrayLike = None, data: jax.typing.ArrayLike = None)

Bases: BaseRenderer

Class to render sources based on rendering them in Fourier space. Sersic profiles are modeled as a series of Gaussian following Shajib (2019) (https://arxiv.org/abs/1906.08263) and the implementation in lenstronomy (https://github.com/lenstronomy/lenstronomy/blob/main/lenstronomy/LensModel/Profiles/gauss_decomposition.py)

frac_start: float

frac_end: float

n_sigma: int

precision: int

use_interp_amps: bool

etas: jax.numpy.array

betas: jax.numpy.array

n_ax: jax.numpy.array

amps_n_ax: jax.numpy.array

get_amps_sigmas(flux, r_eff, n)

render_sersic_mog_fourier(xc, yc, flux, r_eff, n, ellip, theta)

render_sersic(params: dict) → jax.numpy.array

Render a Sersic profile

Parameters:

xc (float) – Central x position
yc (float) – Central y position
flux (float) – Total flux
r_eff (float) – Effective radius
n (float) – Sersic index
ellip (float) – Ellipticity
theta (float) – Position angle in radians

Returns:

Rendered Sersic model

Return type:

jax.numpy.array

render_pointsource(params: dict) → jax.numpy.array

Render a Point source

Parameters:

xc (float) – Central x position
yc (float) – Central y position
flux (float) – Total flux

Returns:

rendered pointsource model

Return type:

jax.numpy.array

Bases: FourierRenderer

Class to render sources based on the hybrid rendering scheme introduced in Lang (2020). This avoids some of the artifacts introduced by rendering sources purely in Fourier space. Sersic profiles are modeled as a series of Gaussian following Shajib (2019) (https://arxiv.org/abs/1906.08263) and the implementation in lenstronomy (https://github.com/lenstronomy/lenstronomy/blob/main/lenstronomy/LensModel/Profiles/gauss_decomposition.py).

Our scheme is implemented slightly differently than Lang (2020), specifically in how it chooses which gaussian components to render in Fourier vs. Real space. Lang (2020) employs a cutoff based on distance to the edge of the image. However given some of jax’s limitation with dynamic shapes (see more here -> https://jax.readthedocs.io/en/latest/notebooks/Common_Gotchas_in_JAX.html#dynamic-shapes), we have not implemented that specific criterion. Instead we use a simpler version where the user must decide how many components to render in real space, starting from the largest ones. While this is not ideal in all circumstances it still overcomes many of the issues of rendering purely in fourier space discussed in Lang (2020).

num_pixel_render: int

w_real: jax.numpy.array

w_fourier: jax.numpy.array

sig_psf_approx: float

render_sersic_hybrid(xc, yc, flux, r_eff, n, ellip, theta)

render_sersic(params: dict) → jax.numpy.array

Render a Sersic profile

Parameters:

xc (float) – Central x position
yc (float) – Central y position
flux (float) – Total flux
r_eff (float) – Effective radius
n (float) – Sersic index
ellip (float) – Ellipticity
theta (float) – Position angle in radians

Returns:

Rendered Sersic model

Return type:

jax.numpy.array

render_pointsource(params: dict) → jax.numpy.array

Render a Point source

Parameters:

xc (float) – Central x position
yc (float) – Central y position
flux (float) – Total flux

Returns:

rendered pointsource model

Return type:

jax.numpy.array

class pysersic.PixelRenderer(im_shape: Iterable, pixel_PSF: jax.numpy.array, os_pixel_size: int | None = 6, num_os: int | None = 12)

Bases: BaseRenderer

Render class based on rendering in pixel space and then convolving with the PSF

os_pixel_size: int

num_os: int

w_os: jax.numpy.array

x_os_lo: int

x_os_hi: int

y_os_lo: int

y_os_hi: int

X_os: jax.numpy.array

Y_os: jax.numpy.array

render_int_sersic(xc, yc, flux, r_eff, n, ellip, theta)

render_sersic(params: dict) → jax.numpy.array

Render a sersic profile

Parameters:

xc (float) – Central x position
yc (float) – Central y position
flux (float) – Total flux
r_eff (float) – Effective radius
n (float) – Sersic index
ellip (float) – Ellipticity
theta (float) – Position angle in radians

Returns:

Rendered Sersic model

Return type:

jax.numpy.array

render_pointsource(params: dict) → jax.numpy.array

Render a Point source by interpolating given PSF into image. Currently jax only supports linear intepolation.

Parameters:

xc (float) – Central x position
yc (float) – Central y position
flux (float) – Total flux

Returns:

rendered pointsource model

Return type:

jax.numpy.array

class pysersic.PySersicResults(data: ArrayLike, rms: ArrayLike, psf: ArrayLike, loss_func: Callable, renderer: pysersic.rendering.BaseRenderer, mask: ArrayLike | None = None)

__repr__() → str: Return repr(self).

add_method_used(method)

add_prior(prior)

add the prior object to the result object once created

Parameters:: prior (pysersic.priors.BasePrior) – created prior object

Method to injest data from optimized SVI model or results of sampling. When sampling data is input, sets the class attribute sampling_results, while for an SVI run, saves svi_results. A PysersicResults object can have at most data for one SVI and one sampling run. Both class atttributes are instances of arviz.InferenceData objects.

Parameters:

sampler (Optional[numpyro.infer.mcmc.MCMC], optional) – numpyro sampler containing results
svi_res_dict (Optional[dict], optional) – Dictionary containing ‘guide’, ‘model’ and ‘svi_result’ specifying a trained SVI model
purge_extra (Optional[bool], optional) – Whether to purge variables containing ‘auto’, ‘base’ or ‘unwrapped’ often used in reparamaterization, by default True
num_sample (Optional[int]) – Number of samples to draw from trained SVI posterior, no effect if sampling was used.
rkey (Optional[jax.random.PRNGKey], optional) – PRNG key to use, by default jax.random.PRNGKey(5)

Returns:

ArviZ Summary of results

Return type:

pandas.DataFrame

Raises:

AssertionError – Must supply one of sampler or svi_dict

_parse_injested_data(data: arviz.InferenceData, purge_extra: bool = True, save_model: bool = True) → arviz.InferenceData

Helper function to postprocess the poterior object (internal use).

Parameters:

data (arviz.InferenceData) – _description_
purge_extra (bool, optional) – whether to purge extra params not part of the fitting, by default True
save_model (bool) – Whether to set self.models with model images from posterior

Returns:

the cleaned up object

Return type:

arviz.InferenceData

summary() → pandas.DataFrame

Convenience function for returning the summary dataframe using the arviz summary.

Returns:: data frame containing the arviz summary of the fit.
Return type:: pandas.DataFrame

get_median_model()

corner(quantiles=[0.16, 0.5, 0.84], **kwargs)

Return a corner plot of the parameter estimation

Parameters:

quantiles (ListLike, optional) – which quantiles to mark on the corner plot, can pass None, by default [0.16,0.5,0.84]
**kwargs – any additional arguments to pass to corner (for the single plot case).

Returns:

fig object containing the corner plots

Return type:

matplotlib.figure

retrieve_param_quantiles(quantiles: ListLike = [0.16, 0.5, 0.84], return_dataframe: bool = False) → pandas.DataFrame | dict

retrieve quantiles on the parameter estimation

Parameters:

quantiles (ListLike, optional) – array of quantiles to pull, must be between 0 and 1, by default [0.16,0.5,0.84]
return_dataframe (bool, optional) – whether to return dataframe instead of simple dict, by default False

Returns:

dict or dataframe with index/keys as parameters and columns/values as the chosen quantiles.

Return type:

Union[pd.DataFrame,dict]

retrieve_med_std(return_dataframe: bool = False) → pandas.DataFrame | dict

latex_table(quantiles: ListLike = [0.16, 0.5, 0.84])

Generate a simple AASTex deluxetable with the fit parameters. Prints the result.

Parameters:: quantiles (ListLike, optional) – quantiles to use must be len 3 as we do upper-median and median-lower to get +/- values, by default [0.16,0.5,0.84]
Raises:: AssertionError – if the quantile list does not have three values

get_chains() → xarray.Dataset

Wrapper for az.extract, producing the chains/draws for the run

Returns:: chain object
Return type:: xarray.Dataset

compute_statistic(parameter: str, func: Callable) → ArrayLike

Compute an arbitrary array statistic on the chain for a given parameter. For example, the std of all ellipticity draws. or the mean of the sersic draws

Parameters:

parameter (str) – a legal parameter from the fit. Use e.g. results.svi_summary() to see them.
func (Callable) – any function which reads in array-like data and computes something

Returns:

the computed statistic(s)

Return type:

ArrayLike

save_result(fname: str)

Save summary of the fit, copies of the data, the chains, and some other info about priors and renderers into an asdf file for later retrieval.

Parameters:: fname (str) – filename to save to.

sample_posterior(num_sample: int, purge_extra: bool | None = True, rkey: jax.random.PRNGKey | None = random.PRNGKey(7)) → arviz.InferenceData

Generate extra samples from an trained SVI posterior

Parameters:

num_sample (int) – number of samples to draw
purge_extra (Optional[bool], optional) – Whether to purge variables containing ‘auto’, ‘base’ or ‘unwrapped’ often used in reparamaterization, by default True
rkey (Optional[jax.random.PRNGKey], optional) – PRNG key to use, by default jax.random.PRNGKey(7)

Returns:

arviz InferenceData object containing posterior

Return type:

az.InferenceData

pysersic.parse_multi_results(results: PySersicResults, source_num: int) → PySersicResults

Function written to parse results from a FitMulti instance to isolate a single source. A new PySersicResults class is created with only the posteriors of the specified source. The original chains saved under .idata_all

Parameters:

results (PySersicResults) – Results class from a FitMulti instance
source_num (int) – Source number to isolate or if equal to -1 will reset to the joint posterior of all sources

Returns:

Results class with all the meta-data the same but the specified source isolated. The original posterior is saved under .idata_all

Return type:

PySersicResults