Inline Analysis

When running a large number of models, each of which takes a few seconds (or more), it’s important to do as much analysis “inline” as possible. For example, say we are interested in obtaining confidence intervals for quantities other than the free parameters of our model. Yes, we could go back later and re-run certain subsets of models and extract whatever information we want, but with a little planning, we can eliminate the need for these “extra” computations. The emcee code dubs such quantities arbitrary meta-data blobs, and as a result, any quantities computed during calculations in ares will be named “blobs” as well.

The goal of this section is to outline the general procedure used to save meta-data blobs of your choosing, which can be tricky because different quantities of interest are computed in very different ways and often are of diverse shapes and variables types.

A few examples of meta-data blobs:

  • Scalar blobs (e.g., the CMB optical depth, \(\tau_e\), the midpoint of reionization, \(z_{\mathrm{rei}}\))
  • 1-D blobs (e.g., the global 21-cm signal, \(\delta T_b(z)\), the thermal history, \(T_K(z)\))
  • 2-D blobs (e.g., the star-formation efficiency, \(f_{\ast}(z, M_h)\), the meta-galactic radiation background intensity, \(J(z, \nu)\))

Example: Scalar and 1-D blobs only

Let’s learn by example. Here is a typical calculation (model for the global 21-cm signal), where we have modified the input parameters so that a few quantities of interest are saved:

  • Scalars: extrema in the global 21-cm signal (which we label turning points B, C, and D).
  • 1-D arrays: the full ionization history, thermal history, and the global 21-cm signal.

So, we define a nested list containing the names of our blobs:

blob_names = [['tau_e'], ['cgm_h_2', 'igm_Tk', 'igm_dTb']]

The blobs ares sorted by their dimensionality: the first sublist contains the names of all scalar blobs, while the second contains the 1-D blobs. Important question: how do you know the names of blobs? The scalar blobs, in this case just tau_e (the CMB optical depth) are all attributes of the Global21cm simulation class (well, really of ares.analysis.Global21cm, but they get inherited). The 1-D blobs are all names of ares fields: see Field Listing for more information.

Now, for the 1-D blobs we also need to provide a sequence of redshifts at which to save each quantity:

blob_ivars = [None, np.arange(6, 21)]

Notice that blob_ivars is a 2-element list (ivars is short for “independent variables,” since in general they need not be redshifts): one element for each blob group (scalar and 1-D). Since the scalars are just numbers, the first element in this list is just None, while the second indicates that we’ll save the desired quantities at redshifts \(z=6,7,...,20\).

We supply these lists via parameters of the same name:

pars = \
{
 'problem_type': 101,           # Simple global 21-cm problem
 'blob_names': blob_names,
 'blob_ivars': blob_ivars,
 'tanh_model': True,            # Just to speed things up
}

Now, we just run the simulation in the usual way:

sim = ares.simulations.Global21cm(**pars)
sim.run()

sim.GlobalSignature()

To verify that the simulation knows about our blobs, some basic information is available via attributes:

print sim.blob_names
print sim.blob_ivars
print sim.blob_dims
print sim.blob_nd

The blobs themselves can be accessed via:

sim.blobs

In this case, using blobs isn’t necessary since we have all data from the simulation at our fingertips. However, the attribute sim.blobs is extremely important for model grids or MCMC fits, as its contents are the only thing written to disk other than the MCMC samples themselves.

Special Redshifts

We’ll often be interested in saving a series of quantities at the redshifts corresponding to the extrema of the global 21-cm signal, or the midpoint of reionization, etc. However, since those aren’t known a-priori, we can’t specify them like we did above. Instead, we tag a suffix (either B, C, or D) onto pre-existing ares fields, i.e.,

extrema = ['z_B', 'igm_dTb_B', 'z_C', 'igm_dTb_C', 'igm_Tk_C']
blob_names = [extrema, ['cgm_h_2', igm_Tk', 'igm_dTb']]

and so on.

Example: 2-D blobs

Now, let’s track slightly more complex blobs. For example, if you’re running models of galaxy populations (see More Detailed Models: The GalaxyPopulation object), you might want to save the galaxy luminosity function at a series of magnitudes and a series of redshifts.