Analyze Model¶

Class to aid in GWB model analysis

class gwb_analysis.analyze_model.Model_Info(path, file, model_name, color, line_style, threshold=0.5, evolving=False, stdev=None, nfreq=5, param_space_name='PS_Astro_Strong_All')¶

Bases: object

Class to analyze gravitational wave background models generated using holodeck. This class automatically reads the data and assigns attributes for the GWB amplitudes, parameters, likelihoods, and frequencies.

See also

This class works with outputs from holodeck.

Parameters:

pathstr: Path to the data, should end with a ‘/’ and be the same as the path used to generate the data
filestr: Filename of the data, should be the same as the filename used to generate the data
model_namestr: Name of the model, used for plotting
colorstr: Color for all plotting associated with this model
line_stylestr: Line style for all plotting associated with this model
thresholdfloat, optional: Minimum likelihood cutoff, used for plotting error bars on the spectrum. Default is 0.5
evolvingbool, optional: Whether the model has MMBulge evolution, used for plotting. Default is False
stdevfloat, optional: Standard deviation of evolution parameter, used for plotting, should be None if evolving is False. Default is None
nfreqint, optional: Number of frequency bins to use in likelihood calculation, 5 is recommended. Default is 5
param_space_namestr, optional: Name of parameter space, should be the same as the parameter space used to generate the data. Default is ‘PS_Astro_Strong_All’

Attributes:

pathstr: Path to the data specified as an argument
filestr: Filename of the data specified as an argument
param_space_namestr: Name of parameter space specified as an argument
gwbarray-like: GWB amplitudes for each model, shape is (nsamples, nrealizations, nfreqs)
param_namesarray-like: Names of parameters in the model
paramsarray-like: Parameter sample associated with each model
ln_likearray-like: Log-likelihood of each model
freqsarray-like: Frequencies at which the amplitudes are calculated
model_namestr: Model name specified as an argument
colorstr: Color for all plotting associated with this model specified as an argument
line_stylestr: Line style for all plotting associated with this model specified as an argument
thresholdfloat: minimum likelihood cutoff specified as an argument
evolvingbool: Whether the model has MMBulge evolution specified as an argument
stdevfloat: Standard deviation of evolution parameter specified as an argument
nfreqint: Number of frequency bins to use in likelihood calculation specified as an argument
idcsarray-like: Index of evolution parameter, None if evolving is False
space_classclass: Class of the parameter space, used to get priors and fiducial values
fiducial_values: dict: Dictionary of fiducial values for each parameter in the model
paramsdict: Dictionary of median posterior values for each parameter in the model, will be the same as fiducial unless get_posteriors() is called
params_errdict: Dictionary of standard deviation of posterior values for each parameter in the model, will be the same as fiducial values unless get_posteriors() is called
plt_labelsdict: Dictionary of plot labels associated with each parameter, used for plotting

Methods

`L_from_Mbh_via_lambda`(mbh_log10, knee, norm, ...)	Calculate AGN luminosity function by convolving black hole mass function with an Eddington ratio distribution function
`L_from_Mbh_via_mdot_eta_func`(mbh_log10, ...)	Calculate luminosity from black hole mass using the accretion rate and radiative efficiency.
`add_spectrum`(ax[, lw, errorbars, label])	Add the spectrum of the model to a given axis, with optional error bars and label.
`bhar_bondi`(mbh_log10[, rho, cs, mth])	Calculate black hole accretion rate from bondi accretion.
`bhar_gal`(mbh_log10, redshift[, mth])	Calculate black hole accretion rate as a function of black hole mass and redshift using the MMBulge relation and the GSMF.
`bhmf`(mbh_log10, redz)	Produce the black hole mass function at a given redshift.
`bhmf_err`(mbh_log10, redz[, ndraws])	Produces the black hole mass function at a given redshift.
`bhmf_from_gsmf`(mstar_log10, mbh_log10, redshift)	Like bhmf_conv in holodeck except this starts with a GSMF and the convolution is done via dot product
`calculate_radiative_efficiency`(zval, mbh_log10)	Calculate the radiative efficiency implied by the model at a given redshift by comparing the change in the black hole mass function between two redshifts to the luminosity function at the average redshift.
`corner_plot`([nbins, cmap])	Old and slow, but tested version of corner plot, will be removed in favor of corner_plot_fast() pending appropriate testing.
`corner_plot_fast`([nbins, cmap])	Faster version of corner plot.
`eta_from_mbh_davis`(mbh_log10)	Calulate radiative efficiency as a function of black hole mass using the fit from Davis & Laor (2011).
`eta_from_mbh_line`(mbh_log10[, mth])	Calculate radiative efficiency as a function of black hole mass using the a line fit to the data from Li et al. (2012).
`eta_from_mbh_logistic`(mbh_log10[, min, k, ...])	Calulate radiative efficiency as a function of black hole mass a logistic fit to the data from Li et al. (2012).
`fdfunc`(mbh_log10, redshift[, fdmin])	Fit to agn fraction as a function of stellar mass from Zou et al. (2024).
`get_posteriors`()	Get the median posterior values for each parameter in the model
`get_priors`()	Returns parameter names and sampled distributions from the prior parameter space
`get_shend`(redshift[, path_to_shen_data])	Retrieve the data from the Shen et al. 2020 bolometric luminosity function at a given redshift.
`get_shenf`(redshift[, path_to_shen_fits])	Retrieve the fit to the Shen et al. 2020 bolometric luminosity function at a given redshift.
`gsmf`(mstar_log10, redz)	Produces the galaxy stellar mass function at a given redshift.
`loglam_func`(mbh_log10, knee, norm, slope[, ...])	Calculate Eddington fraction as a function of mass.
`loglam_func_line`(mbh_log10[, mth])	Calculate Eddington fraction as a function of mass.
`plot_histogram`(ax, param_name[, nbins, ...])	Plot a histogram of the posterior distribution for a given parameter, with optional prior distribution overlaid.

get_priors()¶

Returns parameter names and sampled distributions from the prior parameter space

Sampled distribution of parameter i: space.param_samples.transpose()[i] Name of parameter i: space.param_names[i]

get_posteriors()¶: Get the median posterior values for each parameter in the model

plot_histogram(ax, param_name, nbins=20, histtype='bar', label=None, prior=False)¶

Plot a histogram of the posterior distribution for a given parameter, with optional prior distribution overlaid.

Parameters:

axmatplotlib axis: axis to which the histogram will be added
param_namestr: name of the parameter to plot, should be one of the parameter names in the model
nbinsint, optional: number of bins to use in the histogram, default is 20
labelstr, optional: label for the histogram, default is None, in which case the parameter name will be used
priorbool, optional: whether to plot the prior distribution overlaid on the histogram, default is False

Returns:

None, the histogram is added to the given axis

corner_plot(nbins=20, cmap='Blues')¶

Old and slow, but tested version of corner plot, will be removed in favor of corner_plot_fast() pending appropriate testing. Creates a corner plot for all parameters in the model.

Parameters:

nbinsint, optional: Number of bins to use in the histograms and contour plots. Default is 20
cmapstr, optional: Name of matplotlib colormap to use for the points in the corner plot. Default is ‘Blues’

Returns:

None, the corner plot is displayed using matplotlib

corner_plot_fast(nbins=20, cmap='Blues')¶

Faster version of corner plot. Creates a corner plot for all parameters in the model.

Parameters:

nbinsint, optional: Number of bins to use in the histograms and contour plots. Default is 20
cmapstr, optional: Name of matplotlib colormap to use for the points in the corner plot. Default is ‘Blues’

Returns:

None, the corner plot is displayed using matplotlib

add_spectrum(ax, lw=3, errorbars=False, label=None)¶

Add the spectrum of the model to a given axis, with optional error bars and label.

Parameters:

axMatplotlib axis: axis to which the spectrum will be added
lwfloat, optional: Line width of the spectrum, default is 3
errorbarsbool, optional: Whether to add error bars to the spectrum, default is False
labelstr, optional: Label for the spectrum, default is None, in which case the model name will be used

Returns:

None, the spectrum is added to the given axis

bhmf(mbh_log10, redz)¶

Produce the black hole mass function at a given redshift. This is calculated using holodeck by convolving a double Schechter GSMF with an Mbh–M_bulge relation \(M_\mathrm{BH} = \alpha_0 \left(\frac{M_{\mathrm{bulge}}}{10^{11}\, M_{\odot}}\right)^{\beta_0}\)

Parameters:

massarray: Black hole masses at which the BHMF is evaluated, in log10(Mbh/Msol)
redzfloat: Redshift at which the BHMF is evaluated

Returns:

bhmfarray-like: Black hole mass function at the given redshift

bhmf_err(mbh_log10, redz, ndraws=100)¶

Produces the black hole mass function at a given redshift. Calculated using holodeck by connvolving a double Schechter GSMF with a MMBulge relation

Parameters:

massarray: Black hole masses at which to evaluate the BHMF, in log10(Mbh/Msol)
redzfloat: Redshift at which the BHMF is evaluated

Returns:

bhmf_medianarray: The median black hole mass function at the given redshift, calculated using the median posterior values for the parameters in the model
bhmf_upper_boundarray: The 84th percentile of the black hole mass function at the given redshift
bhmf_lower_boundarray: The 16th percentile of the black hole mass function at the given redshift

gsmf(mstar_log10, redz)¶

Produces the galaxy stellar mass function at a given redshift. Calculated using holodeck by connvolving a double Schechter GSMF with a MMBulge relation

Parameters:

massarray: Black hole masses at which to evaluate the BHMF, in log10(Mbh/Msol)
redzfloat: Redshift at which the BHMF is evaluated

Returns:

gsmfarray: The galaxy stellar mass function at the given redshift

get_shenf(redshift, path_to_shen_fits='/Users/cayenne/Documents/Research/quasarlf/qlffits/')¶

Retrieve the fit to the Shen et al. 2020 bolometric luminosity function at a given redshift.

Parameters:

redshiftflaot: Redshift of the fit to be used 0.2-7.0 in steps of 0.2
path_to_shen_fitsstr, optional: Path to the fits. Default is “/Users/cayenne/Documents/Research/quasarlf/qlffits/”

Returns:

xarray: The x-axis values, luminosity in log10(L/erg/s)
yarray: Fit to log phiL

get_shend(redshift, path_to_shen_data='/Users/cayenne/Documents/Research/quasarlf/qlfdata/')¶

Retrieve the data from the Shen et al. 2020 bolometric luminosity function at a given redshift.

Parameters:

redshiftflaot: Redshift of the fit to be used 0.2-7.0 in steps of 0.2
path_to_shen_fitsstr, optional: Path to the data. Default is “/Users/cayenne/Documents/Research/quasarlf/qlfdata/”

Returns:

xarray: The x-axis values, luminosity in log10(L/erg/s)
yarray: Data for log phiL

calculate_radiative_efficiency(zval, mbh_log10, step=0.001, path_to_shen_fits='/Users/cayenne/Documents/Research/quasarlf/qlffits/')¶

Calculate the radiative efficiency implied by the model at a given redshift by comparing the change in the black hole mass function between two redshifts to the luminosity function at the average redshift. This is a rough calculation that assumes that the change in the BHMF and LF between the two redshifts is solely due to accretion and does not consider mergers or other processes that may contribute to the growth of black holes.

May have bugs, shouldn’t be used

Parameters:

zvalfloat: The redshift at which to calculate the radiative efficiency
mbh_log10array-like: The log10 of the black hole masses at which to evaluate the BHMF.
stepfloat, optional: The step in redshift to use for calculating the change in the BHMF and LF. Default is 1e-3
path_to_shen_fits: str, optional: Path to the Shen et al. 2020 fits for the bolometric LF at different redshifts. Default is “/Users/cayenne/Documents/Research/quasarlf/qlffits/”

Returns:

eradfloat: Radiative efficiency between the two redshifts
mdotfloat: The total integrated mass density gain between those redshifts (scaled)
Lumfloat: The total integrated luminosity at the latter redshift

Warning

Does not incorporate AGN Fraction
Radiative efficiency is a constant here
Need to make sure that it is always integrating over roughly similar mass - luminosity ranges
Integration only cosiders two bins and nothing in between, but should be fine for order of magnitude calculation
May have a volume normalization issue

fdfunc(mbh_log10, redshift, fdmin=0.0)¶

Fit to agn fraction as a function of stellar mass from Zou et al. (2024). Here stellar mass is inferred from black hole mass

Parameters:

mbh_log10array-like: Log10 of black hole mass in solar masses
redshiftfloat: Redshift at which to evaluate the AGN fraction
fdminfloat, optional: Minimum AGN fraction to return, default is 0.03, which is the value used in Shen et al. 2020

Returns:

phi_fdarray-like: AGN fraction as a function of black hole mass at the given redshift

bhmf_from_gsmf(mstar_log10, mbh_log10, redshift)¶

Like bhmf_conv in holodeck except this starts with a GSMF and the convolution is done via dot product

Parameters:

mstar_log10array-like: Log10 of stellar mass in solar masses at which to evaluate the BHMF
mbh_log10array-like: Log10 of black hole mass in solar masses at which to evaluate the BHMF
redshiftfloat: Redshift at which to evaluate the BHMF

Returns:

bhmf_convarray-like: The black hole mass function at the given redshift calculated from the GSMF and MMBulge relation

bhar_gal(mbh_log10, redshift, mth=None)¶

Calculate black hole accretion rate as a function of black hole mass and redshift using the MMBulge relation and the GSMF.

The function in the paper is in terms of stellar mass, but we can use the MMBulge relation to convert it to a function of black hole mass. Fit from Zou et al. (2024).

Msun / year Error ranges from 0.1 - 0.3 dex depending on redshift and mass (see Figure 6)

Parameters:

mbh_log10array-like: Log10 of black hole mass in solar masses at which to evaluate the BHAR
redshiftfloat: Redshift at which to evaluate the BHAR
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

bhararray-like: Black hole accretion rate in Msun / year as a function of black hole mass and redshift

bhar_bondi(mbh_log10, rho=0.1, cs=100, mth=None)¶

Calculate black hole accretion rate from bondi accretion.

Msun / year Error ranges from 0.1 - 0.3 dex depending on redshift and mass (see Figure 6)

Parameters:

mbh_log10array-like: Log10 of black hole mass in solar masses at which to evaluate the BHAR
rhofloat: Density of the gas in the vicinity of the black hole, in units of g / cm^3, Default is 0.1
csfloat: Sound speed of the gas in the vicinity of the black hole, in units of km/s, Default is 100
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

bhararray-like: Black hole accretion rate in Msun / year as a function of black hole mass and redshift

eta_from_mbh_davis(mbh_log10)¶

Calulate radiative efficiency as a function of black hole mass using the fit from Davis & Laor (2011).

Parameters:

mbh_log10array-like: Log10 of black hole mass in solar masses
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

etasarray-like: Radiative efficiency as a function of black hole mass

eta_from_mbh_line(mbh_log10, mth=None)¶

Calculate radiative efficiency as a function of black hole mass using the a line fit to the data from Li et al. (2012). Two lines of different slopes with a cutoff.

Caution

Not advised to use for redshifts below 0.8.

Parameters:

mbh_log10array-like: Log10 of black hole mass in solar masses at which to evaluate the radiative efficiency
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

etasarray-like: Radiative efficiency as a function of black hole mass

eta_from_mbh_logistic(mbh_log10, min=0.001, k=-1.4, m0=9.7, mth=None)¶

Calulate radiative efficiency as a function of black hole mass a logistic fit to the data from Li et al. (2012). No redshift dependence.

Tip

Use m0 = 8.4 for a smoother transition between high and low values of radiative efficiency.

Caution

Not advised to use for redshifts below 0.8.

Parameters:

mbh_log10array-like: Log10 of black hole mass in solar masses at which to evaluate the radiative efficiency
minfloat, optional: Minimum value of the logistic function, default is 0.001. Default is 0.001
kfloat, optional: Stepth of the logistic function. Default is -1.4
m0float, optional: The value of mbh_log10 at which the logistic function is halfway between its minimum and maximum values. Default is 9.1
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

etasarray-like: Radiative efficiency as a function of black hole mass

L_from_Mbh_via_mdot_eta_func(mbh_log10, lums_log10, redshift, ndens=None, scatter=None, eta_func='Davis', rad_eff=None, mdot_func='Gal', mth=None)¶

Calculate luminosity from black hole mass using the accretion rate and radiative efficiency. The accretion rate is calculated using the MMBulge relation and the GSMF, and the radiative efficiency is calculated using one of several functions of black hole mass.

Parameters:

mbh_log10array: Black hole masses to evaluate the luminosity function at, in log10(Mbh/Msol)
lums_log10array: Array of log10 luminosities at which to evaluate the luminosity function
redshiftfloat: Redshift at which to evaluate the luminosity function
ndensarray, optional: Number density of black holes at the given masses and redshift. If not provided, it will be calculated using the bhmf function. Default is None.
eta_funcbool, optional: Which functional form to use for calculating radiative efficiency, options are ‘Davis’, ‘Logistic’, ‘Line’, and ‘Constant’. Default is ‘Davis’
rad_efffloat, optional: The constantvalue of the radiative efficiency to use when eta_func is ‘Constant’. Default is None
mdot_funcbool, optional: Which functional form to use for calculating accretion rate, options are ‘Gal’, ‘Bondi’, and ‘Lambda’. Default is ‘Gal’.
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

lf_convarray: The luminosity function calculated from the black hole mass function, accretion rate, and radiative efficiency

Raises:

ValueError: If eta_func is ‘Constant’ and rad_eff is not provided, a ValueError is raised.

loglam_func(mbh_log10, knee, norm, slope, lowlam=-5, hilam=1, mth=None)¶

Calculate Eddington fraction as a function of mass. Schechter function

Parameters:

mbh_log10array: Black hole masses to evaluate the luminosity function at, in log10(Mbh/Msol)
kneefloat: The turnover point of the Schechter function
normfloat: The normalization of the Schechter function
slopefloat: The slope of the Schechter function
lowlamfloat: Lower limit of allowable Eddington ratios
hilamfloat: Upper limit of allowable Eddington ratios
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

loglam_Marray: Log10 of the Eddington fraction as a function of black hole mass. Functional form is a Schechter function with the given knee, norm, and slope, and is clipped to be between lowlam and hilam for numerical reasons.

loglam_func_line(mbh_log10, mth=None)¶

Calculate Eddington fraction as a function of mass. Schechter function

Parameters:

mbh_log10array: Black hole masses to evaluate the luminosity function at, in log10(Mbh/Msol)
kneefloat: The turnover point of the Schechter function
normfloat: The normalization of the Schechter function
slopefloat: The slope of the Schechter function
lowlamfloat: Lower limit of allowable Eddington ratios
hilamfloat: Upper limit of allowable Eddington ratios
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

loglam_Marray: Log10 of the Eddington fraction as a function of black hole mass. Functional form is a Schechter function with the given knee, norm, and slope, and is clipped to be between lowlam and hilam for numerical reasons.

L_from_Mbh_via_lambda(mbh_log10, knee, norm, slope, sigma_loglam, redshift, logL_grid, ndens=None, loglam_func='Schechter', lowlam=-15, hilam=11, mth=None)¶

Calculate AGN luminosity function by convolving black hole mass function with an Eddington ratio distribution function

Parameters:

mbh_log10array: Black hole masses to evaluate the luminosity function at, in log10(Mbh/Msol)
kneefloat: The turnover point of the Schechter function
normfloat: The normalization of the Schechter function
slopefloat: The slope of the Schechter function
sigma_loglamfloat: The scatter in log lambda at fixed black hole mass, which is assumed to be Gaussian
redshiftfloat: Redshift at which to evaluate the luminosity function
logL_gridarray: The grid of log luminosities at which to evaluate the luminosity function
ndensarray, optional: Number density of black holes at the given masses and redshift. If not provided, it will be calculated using the bhmf function. Default is None.
loglam_funcstr, optional: Which functional form to use for calculating the mean log lambda as a function of black hole mass, options are ‘Schechter’ and ‘Line’. Default is ‘Schechter’
lowlamfloat: Lower limit of allowable Eddington ratios
hilamfloat: Upper limit of allowable Eddington ratios
mthmodule, optional: Module to use for mathematical functions, default is numpy.

Returns:

lum_funcarray: The luminosity function calculated from the black hole mass function and Eddington ratio distribution function