# PIG 21 - Models improvements#

Author: Axel Donath, Régis Terrier and Quentin Remy

Created: Jun 10, 2020

Accepted: July 31, 2020

Status: accepted

Discussion: GH 2944

## Abstract#

This PIG outlines further minor improvements to the modeling framework of Gammapy. This includes introduction of a spectra models with norm parameters, minimal support for energy dependent spatial models as well as simplifications of the YAML serialisation.

## Proposal#

### Spectral Norm Models#

For the purpose of adjusting template based models we propose to introduce a new
class of spectral models. Those spectral models feature a norm parameter instead
of amplitude and are named using the `NormSpectralModel`

suffix:

```
from gammapy.modeling.models import (
PowerLawNormSpectralModel,
LogParabolaNormSpectralModel,
PiecewiseBrokenPowerlawNormSpectralModel,
ConstantNormSpectralModel,
)
pwl_norm = PowerLawNormSpectralModel(norm=, tilt=, reference=)
log_parabola_norm = LogParabolaNormSpectralModel(norm=, alpha=, beta=)
bpwl_norm = PiecewiseBrokenPowerlawNormSpectralModel(norms=, energies=)
const_norm = ConstantNormSpectralModel(norm=)
# typically used like
template = TemplateSpectralModel()
model = template * pwl_norm
model = BackgroundModel(
map=map,
spectral_model=pwl_norm
)
model = SkyDiffuseCube(
map=map,
spectral_model=pwl_norm
)
```

This requires removing the hard-coded parameters of the `TemplateSpectralModel`

and `BackgroundModel`

.

### Energy Dependent Spatial Models#

A very common use-case in scientific analyses is to look for energy dependent morphology of extended sources. In the past this kind of analysis has been typically done by splitting the data into energy bands and doing individual fits of the morphology in these energy bands. In a combined spectral and spatial (“cube”) analysis this can be naturally achieved by allowing for an energy dependent spatial model, where the energy dependence is e.g. described by a parametric model expression with energy dependent parameters.

In the current model scheme we use a “factorised” representation of the source model, where the spatial, energy and time dependence are independent model components and not correlated:

To represent energy dependent morphology we propose to introduce energy dependent spatial models:

In general the energy dependence is optional. If the spatial model does not declare an energy dependence it assumes the same morphology for all energies. This also ensures backwards compatibility with the current behaviour.

To limit the implementation effort in this PIG we propose to only add an example energy dependent
custom model to our documentation. We do not propose to introduce general dependence of arbitrary
model parameters for any spatial model, such as `GaussianSpatialModel`

or `DiskSpatialModel`

.
An example of how this can be achieved with a custom implemented model is given below:

```
from gammapy.modeling.models import SpatialModel
from astropy.coordinates.angle_utilities import angular_separation
class MyCustomGaussianModel(SpatialModel):
"""My custom gaussian model.
Parameters
----------
lon_0, lat_0 : `~astropy.coordinates.Angle`
Center position
sigma_1TeV : `~astropy.coordinates.Angle`
Width of the Gaussian at 1 TeV
sigma_10TeV : `~astropy.coordinates.Angle`
Width of the Gaussian at 10 TeV
"""
tag = "MyCustomGaussianModel"
lon_0 = Parameter("lon_0", "0 deg")
lat_0 = Parameter("lat_0", "0 deg", min=-90, max=90)
sigma_1TeV = Parameter("sigma_1TeV", "1 deg", min=0)
sigma_10TeV = Parameter("sigma_10TeV", "0.5 deg", min=0)
@staticmethod:
def evaluate(lon, lat, energy, lon_0, lat_0, sigma_1TeV, sigma_10TeV):
"""Evaluate custom Gaussian model"""
sigmas = u.Quantity([sigma_1TeV, sigma_10TeV])
energy_nodes = [1, 10] * u.TeV
sigma = np.interp(energy, energy_nodes, sigmas)
sep = angular_separation(lon, lat, lon_0, lat_0)
exponent = -0.5 * (sep / sigma) ** 2
norm = 1 / (2 * np.pi * sigma ** 2)
return norm * np.exp(exponent)
@property
def evaluation_radius(self):
"""Evaluation radius (`~astropy.coordinates.Angle`)."""
return 5 * self.sigma_1TeV.quantity
```

### Spectral Absorption Model#

In the current handling of absorbed spectral models we have a very special
`Absorption`

model, which is not a spectral model. To resolve this special
case, we propose to refactor the existing code and handle the absorbed case
using a `CompoundSpectralModel`

. The new implementation is used as follows:

```
from gammapy.modeling.models import EBLAbsorptionNormSpectralModel, PowerLawSpectralModel
absorption = EBLAbsorptionNormSpectralModel.from_reference(
redshift=0.1, alpha_norm=1, reference="dominguez"
)
pwl = PowerLawSpectralModel()
spectral_model = absorption * pwl
assert isinstance(spectral_model, CompoundSpectralModel)
```

In addition we propose to rename `.table_model`

to `.to_template_spectral_model(redshift, alpha_norm)`

.

### Additional Models#

In addition we propose to implement the following models in `gammapy.modeling.models`

:

`SersicSpatialModel`

following the parametrisation of the Astropy`Sersic2D`

model.`BrokenPowerLaw`

with the following parametrisation:

We propose to re-introduce the

`PhaseCurveModel`

to compute mean fluxes over time

### Simplify YAML Representation#

#### Simplify YAML Parameter Representation#

To further simplify the structure of the YAML file we propose to remove parameter properties if they are equivalent to the default values. The current representation looks like:

```
spectral:
type: PowerLawSpectralModel
parameters:
- name: index
value: 2.0
unit: ''
min: .nan
max: .nan
frozen: false
error: 0
- name: amplitude
value: 1.0e-12
unit: cm-2 s-1 TeV-1
min: .nan
max: .nan
frozen: false
error: 0
- name: reference
value: 1.0
unit: TeV
min: .nan
max: .nan
frozen: true
error: 0
```

After the simplification:

```
spectral:
type: PowerLawSpectralModel
parameters:
- name: index
value: 2.0
- name: amplitude
value: 1.0e-12
unit: cm-2 s-1 TeV-1
- name: reference
value: 1.0
unit: TeV
frozen: true
```

## Decision#

This PIG received feedback in form of comments to the PR GH 2944 and via private
Slack messages to the authors. Most concerns were raised against the introduction
of `NormSpectralModels`

. In a following discussion among the lead developers
it was decided to stick with the proposal, but keep the abstraction of `SkyDiffuseCube`

to avoid the case, where a spatial model carries the flux information.
A final review announced on the Gammapy and CC mailing list did not yield
any additional comments. Therefore the PIG was accepted on July 31, 2020.