Skip to content

Model

ELQModel module.

This module provides a class definition for the main functionalities of the codebase, providing the interface with the openMCMC repo and defining some plotting wrappers.

ELQModel dataclass

Class for setting up, running, and post-processing the full ELQModel analysis.

Attributes:

Name Type Description
form dict

dictionary detailing the form of the predictor for the concentration data. For details of the required specification, see parameter.LinearCombinationWithTransform() in the openMCMC repo.
transform dict

dictionary detailing transformations applied to the model components. For details of the required specification, see parameter.LinearCombinationWithTransform() in the openMCMC repo.
model Model

full model specification for the analysis, constructed in self.to_mcmc().
mcmc MCMC

MCMC object containing model and sampler specification for the problem. Constructed from the other components in self.to_mcmc().
n_iter int

number of MCMC iterations to be run.
n_thin int

number of iterations to thin by.
fitted_values ndarray

samples of fitted values (i.e. model predictions for the data) generated during the MCMC sampler. Attached in self.from_mcmc().
Source code in src/pyelq/model.py 
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
@dataclass
class ELQModel:
    """Class for setting up, running, and post-processing the full ELQModel analysis.

    Attributes:
        form (dict): dictionary detailing the form of the predictor for the concentration data. For details of the
            required specification, see parameter.LinearCombinationWithTransform() in the openMCMC repo.
        transform (dict): dictionary detailing transformations applied to the model components. For details of the
            required specification, see parameter.LinearCombinationWithTransform() in the openMCMC repo.
        model (Model): full model specification for the analysis, constructed in self.to_mcmc().
        mcmc (MCMC): MCMC object containing model and sampler specification for the problem. Constructed from the
            other components in self.to_mcmc().
        n_iter (int): number of MCMC iterations to be run.
        n_thin (int): number of iterations to thin by.
        fitted_values (np.ndarray): samples of fitted values (i.e. model predictions for the data) generated during the
            MCMC sampler. Attached in self.from_mcmc().

    """

    form: dict = field(init=False)
    transform: dict = field(init=False)
    model: Model = field(init=False)
    mcmc: MCMC = field(init=False)
    n_iter: int = 1000
    n_thin: int = 1
    fitted_values: np.ndarray = field(init=False)

    def __init__(
        self,
        sensor_object: SensorGroup,
        meteorology: Union[Meteorology, MeteorologyGroup],
        gas_species: GasSpecies,
        background: Background = SpatioTemporalBackground(),
        source_model: SourceModel = Normal(),
        error_model: ErrorModel = BySensor(),
        offset_model: PerSensor = None,
    ):
        """Initialise the ELQModel model.

        Model form is as follows:
        y = A*s + b + d + e
        where:
        - y is the vector of observed concentration data (extracted from the sensor object).
        - A*s is the source contribution (from the source model and dispersion model).
        - b is from the background model.
        - d is from the offset model.
        - e is residual error term and var(e) comes from the error precision model.

        Args:
            sensor_object (SensorGroup): sensor data.
            meteorology (Union[Meteorology, MeteorologyGroup]): meteorology data.
            gas_species (GasSpecies): gas species object.
            background (Background): background model specification. Defaults to SpatioTemporalBackground().
            source_model (SourceModel): source model specification. Defaults to Normal().
            error_model (Precision): measurement precision model specification. Defaults to BySensor().
            offset_model (PerSensor): offset model specification. Defaults to None.

        """
        self.sensor_object = sensor_object
        self.meteorology = meteorology
        self.gas_species = gas_species
        self.components = {
            "background": background,
            "source": source_model,
            "error_model": error_model,
            "offset": offset_model,
        }
        if error_model is None:
            self.components["error_model"] = BySensor()
            warnings.warn("None is not an allowed type for error_model: resetting to default BySensor model.")
        for key in list(self.components.keys()):
            if self.components[key] is None:
                self.components.pop(key)

    def initialise(self):
        """Take data inputs and extract relevant properties."""
        self.form = {}
        self.transform = {}
        component_keys = list(self.components.keys())
        if "background" in component_keys:
            self.form["bg"] = "B_bg"
            self.transform["bg"] = False
        if "source" in component_keys:
            self.transform["s"] = False
            self.form["s"] = "A"
        if "offset" in component_keys:
            self.form["d"] = "B_d"
            self.transform["d"] = False
        for key in component_keys:
            self.components[key].initialise(self.sensor_object, self.meteorology, self.gas_species)

    def to_mcmc(self):
        """Convert the ELQModel specification into an MCMC solver object that can be run.

        Executing the following steps:
            - Initialise the model object with the data likelihood (response distribution for y), and add all the
                associated prior distributions, as specified by the model components.
            - Initialise the state dictionary with the observed sensor data, and add parameters associated with all
                the associated prior distributions, as specified by the model components.
            - Initialise the MCMC sampler objects associated with each of the model components.
            - Create the MCMC solver object, using all of the above information.

        """
        response_precision = self.components["error_model"].precision_parameter
        model = [
            location_scale.Normal(
                "y",
                mean=parameter.LinearCombinationWithTransform(self.form, self.transform),
                precision=response_precision,
            )
        ]

        initial_state = {"y": self.sensor_object.concentration}

        for component in self.components.values():
            model = component.make_model(model)
            initial_state = component.make_state(initial_state)

        self.model = Model(model, response={"y": "mean"})

        sampler_list = []
        for component in self.components.values():
            sampler_list = component.make_sampler(self.model, sampler_list)

        self.mcmc = MCMC(initial_state, sampler_list, self.model, n_burn=0, n_iter=self.n_iter, n_thin=self.n_thin)

    def run_mcmc(self):
        """Run the mcmc function."""
        self.mcmc.run_mcmc()

    def from_mcmc(self):
        """Extract information from MCMC solver class once its has run.

        Performs two operations:
            - For each of the components of the model: extracts the related sampled parameter values and attaches these
                to the component class.
            - For all keys in the mcmc.store dictionary: extracts the sampled parameter values from self.mcmc.store and
                puts them into the equivalent fields in the state

        """
        state = self.mcmc.state
        for component in self.components.values():
            component.from_mcmc(self.mcmc.store)
        for key in self.mcmc.store:
            state[key] = self.mcmc.store[key]

    def plot_log_posterior(self, burn_in_value: int, plot: Plot = Plot()) -> Plot():
        """Plots the trace of the log posterior over the iterations of the MCMC.

        Args:
            burn_in_value (int): Burn in value to show in plot.
            plot (Plot, optional): Plot object to which this figure will be added in the figure dictionary

        Returns:
            plot (Plot): Plot object to which this figure is added in the figure dictionary with
                key 'log_posterior_plot'

        """
        plot.plot_single_trace(object_to_plot=self.mcmc, burn_in=burn_in_value)
        return plot

    def plot_fitted_values(self, plot: Plot = Plot()) -> Plot:
        """Plot the fitted values from the mcmc object against time, also shows the estimated background when possible.

        Based on the inputs it plots the results of the mcmc analysis, being the fitted values of the concentration
        measurements together with the 10th and 90th quantile lines to show the goodness of fit of the estimates.

        Args:
            plot (Plot, optional): Plot object to which this figure will be added in the figure dictionary

        Returns:
            plot (Plot): Plot object to which this figure is added in the figure dictionary with key 'fitted_values'

        """
        plot.plot_fitted_values_per_sensor(
            mcmc_object=self.mcmc, sensor_object=self.sensor_object, background_model=self.components["background"]
        )
        return plot

__init__(sensor_object, meteorology, gas_species, background=SpatioTemporalBackground(), source_model=Normal(), error_model=BySensor(), offset_model=None)

Initialise the ELQModel model.

Model form is as follows: y = As + b + d + e where: - y is the vector of observed concentration data (extracted from the sensor object). - As is the source contribution (from the source model and dispersion model). - b is from the background model. - d is from the offset model. - e is residual error term and var(e) comes from the error precision model.

Parameters:

Name Type Description Default
sensor_object SensorGroup

sensor data.

 required
meteorology Union[Meteorology, MeteorologyGroup]

meteorology data.

 required
gas_species GasSpecies

gas species object.

 required
background Background

background model specification. Defaults to SpatioTemporalBackground().

 SpatioTemporalBackground()
source_model SourceModel

source model specification. Defaults to Normal().

 Normal()
error_model Precision

measurement precision model specification. Defaults to BySensor().

 BySensor()
offset_model PerSensor

offset model specification. Defaults to None.

 None
Source code in src/pyelq/model.py 
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
def __init__(
    self,
    sensor_object: SensorGroup,
    meteorology: Union[Meteorology, MeteorologyGroup],
    gas_species: GasSpecies,
    background: Background = SpatioTemporalBackground(),
    source_model: SourceModel = Normal(),
    error_model: ErrorModel = BySensor(),
    offset_model: PerSensor = None,
):
    """Initialise the ELQModel model.

    Model form is as follows:
    y = A*s + b + d + e
    where:
    - y is the vector of observed concentration data (extracted from the sensor object).
    - A*s is the source contribution (from the source model and dispersion model).
    - b is from the background model.
    - d is from the offset model.
    - e is residual error term and var(e) comes from the error precision model.

    Args:
        sensor_object (SensorGroup): sensor data.
        meteorology (Union[Meteorology, MeteorologyGroup]): meteorology data.
        gas_species (GasSpecies): gas species object.
        background (Background): background model specification. Defaults to SpatioTemporalBackground().
        source_model (SourceModel): source model specification. Defaults to Normal().
        error_model (Precision): measurement precision model specification. Defaults to BySensor().
        offset_model (PerSensor): offset model specification. Defaults to None.

    """
    self.sensor_object = sensor_object
    self.meteorology = meteorology
    self.gas_species = gas_species
    self.components = {
        "background": background,
        "source": source_model,
        "error_model": error_model,
        "offset": offset_model,
    }
    if error_model is None:
        self.components["error_model"] = BySensor()
        warnings.warn("None is not an allowed type for error_model: resetting to default BySensor model.")
    for key in list(self.components.keys()):
        if self.components[key] is None:
            self.components.pop(key)

initialise()

Take data inputs and extract relevant properties.

Source code in src/pyelq/model.py 
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
def initialise(self):
    """Take data inputs and extract relevant properties."""
    self.form = {}
    self.transform = {}
    component_keys = list(self.components.keys())
    if "background" in component_keys:
        self.form["bg"] = "B_bg"
        self.transform["bg"] = False
    if "source" in component_keys:
        self.transform["s"] = False
        self.form["s"] = "A"
    if "offset" in component_keys:
        self.form["d"] = "B_d"
        self.transform["d"] = False
    for key in component_keys:
        self.components[key].initialise(self.sensor_object, self.meteorology, self.gas_species)

to_mcmc()

Convert the ELQModel specification into an MCMC solver object that can be run.

Executing the following steps
  • Initialise the model object with the data likelihood (response distribution for y), and add all the associated prior distributions, as specified by the model components.
  • Initialise the state dictionary with the observed sensor data, and add parameters associated with all the associated prior distributions, as specified by the model components.
  • Initialise the MCMC sampler objects associated with each of the model components.
  • Create the MCMC solver object, using all of the above information.
Source code in src/pyelq/model.py 
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
def to_mcmc(self):
    """Convert the ELQModel specification into an MCMC solver object that can be run.

    Executing the following steps:
        - Initialise the model object with the data likelihood (response distribution for y), and add all the
            associated prior distributions, as specified by the model components.
        - Initialise the state dictionary with the observed sensor data, and add parameters associated with all
            the associated prior distributions, as specified by the model components.
        - Initialise the MCMC sampler objects associated with each of the model components.
        - Create the MCMC solver object, using all of the above information.

    """
    response_precision = self.components["error_model"].precision_parameter
    model = [
        location_scale.Normal(
            "y",
            mean=parameter.LinearCombinationWithTransform(self.form, self.transform),
            precision=response_precision,
        )
    ]

    initial_state = {"y": self.sensor_object.concentration}

    for component in self.components.values():
        model = component.make_model(model)
        initial_state = component.make_state(initial_state)

    self.model = Model(model, response={"y": "mean"})

    sampler_list = []
    for component in self.components.values():
        sampler_list = component.make_sampler(self.model, sampler_list)

    self.mcmc = MCMC(initial_state, sampler_list, self.model, n_burn=0, n_iter=self.n_iter, n_thin=self.n_thin)

run_mcmc()

Run the mcmc function.

Source code in src/pyelq/model.py 
158
159
160
def run_mcmc(self):
    """Run the mcmc function."""
    self.mcmc.run_mcmc()

from_mcmc()

Extract information from MCMC solver class once its has run.

Performs two operations
  • For each of the components of the model: extracts the related sampled parameter values and attaches these to the component class.
  • For all keys in the mcmc.store dictionary: extracts the sampled parameter values from self.mcmc.store and puts them into the equivalent fields in the state
Source code in src/pyelq/model.py 
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
def from_mcmc(self):
    """Extract information from MCMC solver class once its has run.

    Performs two operations:
        - For each of the components of the model: extracts the related sampled parameter values and attaches these
            to the component class.
        - For all keys in the mcmc.store dictionary: extracts the sampled parameter values from self.mcmc.store and
            puts them into the equivalent fields in the state

    """
    state = self.mcmc.state
    for component in self.components.values():
        component.from_mcmc(self.mcmc.store)
    for key in self.mcmc.store:
        state[key] = self.mcmc.store[key]

plot_log_posterior(burn_in_value, plot=Plot())

Plots the trace of the log posterior over the iterations of the MCMC.

Parameters:

Name Type Description Default
burn_in_value int

Burn in value to show in plot.

 required
plot Plot

Plot object to which this figure will be added in the figure dictionary

 Plot()

Returns:

Name Type Description
plot Plot

Plot object to which this figure is added in the figure dictionary with key 'log_posterior_plot'
Source code in src/pyelq/model.py 
178
179
180
181
182
183
184
185
186
187
188
189
190
191
def plot_log_posterior(self, burn_in_value: int, plot: Plot = Plot()) -> Plot():
    """Plots the trace of the log posterior over the iterations of the MCMC.

    Args:
        burn_in_value (int): Burn in value to show in plot.
        plot (Plot, optional): Plot object to which this figure will be added in the figure dictionary

    Returns:
        plot (Plot): Plot object to which this figure is added in the figure dictionary with
            key 'log_posterior_plot'

    """
    plot.plot_single_trace(object_to_plot=self.mcmc, burn_in=burn_in_value)
    return plot

plot_fitted_values(plot=Plot())

Plot the fitted values from the mcmc object against time, also shows the estimated background when possible.

Based on the inputs it plots the results of the mcmc analysis, being the fitted values of the concentration measurements together with the 10th and 90th quantile lines to show the goodness of fit of the estimates.

Parameters:

Name Type Description Default
plot Plot

Plot object to which this figure will be added in the figure dictionary

 Plot()

Returns:

Name Type Description
plot Plot

Plot object to which this figure is added in the figure dictionary with key 'fitted_values'
Source code in src/pyelq/model.py 
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
def plot_fitted_values(self, plot: Plot = Plot()) -> Plot:
    """Plot the fitted values from the mcmc object against time, also shows the estimated background when possible.

    Based on the inputs it plots the results of the mcmc analysis, being the fitted values of the concentration
    measurements together with the 10th and 90th quantile lines to show the goodness of fit of the estimates.

    Args:
        plot (Plot, optional): Plot object to which this figure will be added in the figure dictionary

    Returns:
        plot (Plot): Plot object to which this figure is added in the figure dictionary with key 'fitted_values'

    """
    plot.plot_fitted_values_per_sensor(
        mcmc_object=self.mcmc, sensor_object=self.sensor_object, background_model=self.components["background"]
    )
    return plot