Skip to content

Noise Models

Physics-aware observation-noise profiles for source-dependent, replicate-aware, and feature-dependent GPR uncertainty.

matgpr.noise_models

ObservationNoiseProfile dataclass

Per-observation noise standard deviations for GPR workflows.

Gaussian Process Regression libraries commonly accept observation variance through an alpha vector. This profile stores standard deviations because they are easier to reason about in physical target units, then exposes the corresponding variances through :attr:alpha.

Source code in matgpr/noise_models.py
27
28
29
30
31
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
@dataclass(frozen=True)
class ObservationNoiseProfile:
    """Per-observation noise standard deviations for GPR workflows.

    Gaussian Process Regression libraries commonly accept observation variance
    through an ``alpha`` vector. This profile stores standard deviations because
    they are easier to reason about in physical target units, then exposes the
    corresponding variances through :attr:`alpha`.
    """

    noise_std: Sequence[float] | np.ndarray
    labels: Sequence[str] | np.ndarray | None = None
    component_names: Sequence[str] | np.ndarray | None = None

    def __post_init__(self) -> None:
        noise = to_1d_finite(self.noise_std, "noise_std")
        if np.any(noise < 0):
            raise ValueError("noise_std must be non-negative")
        labels = label_array(self.labels, noise.shape[0], default="observation_noise")
        components = label_array(
            self.component_names,
            noise.shape[0],
            default="noise",
        )

        object.__setattr__(self, "noise_std", noise)
        object.__setattr__(self, "labels", labels)
        object.__setattr__(self, "component_names", components)

    @property
    def n_observations(self) -> int:
        """Number of observations in the profile."""
        return int(self.noise_std.shape[0])

    @property
    def alpha(self) -> np.ndarray:
        """Observation variances suitable for scikit-learn GPR ``alpha``."""
        return self.noise_std**2

    @property
    def variance(self) -> np.ndarray:
        """Alias for observation variance."""
        return self.alpha

    def to_frame(self) -> pd.DataFrame:
        """Return a compact dataframe for reports and debugging."""
        return pd.DataFrame(
            {
                "noise_std": self.noise_std,
                "variance": self.variance,
                "label": self.labels,
                "component": self.component_names,
            }
        )

    def summary(self) -> pd.DataFrame:
        """Summarize noise levels by label."""
        frame = self.to_frame()
        return (
            frame.groupby("label", dropna=False)["noise_std"]
            .agg(["count", "mean", "median", "min", "max"])
            .reset_index()
        )

n_observations property

Number of observations in the profile.

alpha property

Observation variances suitable for scikit-learn GPR alpha.

variance property

Alias for observation variance.

to_frame()

Return a compact dataframe for reports and debugging.

Source code in matgpr/noise_models.py
71
72
73
74
75
76
77
78
79
80
def to_frame(self) -> pd.DataFrame:
    """Return a compact dataframe for reports and debugging."""
    return pd.DataFrame(
        {
            "noise_std": self.noise_std,
            "variance": self.variance,
            "label": self.labels,
            "component": self.component_names,
        }
    )

summary()

Summarize noise levels by label.

Source code in matgpr/noise_models.py
82
83
84
85
86
87
88
89
def summary(self) -> pd.DataFrame:
    """Summarize noise levels by label."""
    frame = self.to_frame()
    return (
        frame.groupby("label", dropna=False)["noise_std"]
        .agg(["count", "mean", "median", "min", "max"])
        .reset_index()
    )

SourceNoiseModel dataclass

Assign observation noise by data source, instrument, method, or paper.

Use this when rows come from sources with different expected reliability, such as experimental versus simulated data, different characterization techniques, or multiple literature papers.

Source code in matgpr/noise_models.py
 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
@dataclass(frozen=True)
class SourceNoiseModel:
    """Assign observation noise by data source, instrument, method, or paper.

    Use this when rows come from sources with different expected reliability,
    such as experimental versus simulated data, different characterization
    techniques, or multiple literature papers.
    """

    source_noise_std: Mapping[object, float]
    default_noise_std: float | None = None
    unknown: str = "error"
    label_prefix: str = "source"

    def profile(self, sources) -> ObservationNoiseProfile:
        """Return per-row noise for source labels."""
        source_values = _to_1d_object(sources, "sources")
        normalized_unknown = self.unknown.lower()
        if normalized_unknown not in {"error", "default"}:
            raise ValueError("unknown must be 'error' or 'default'")
        if normalized_unknown == "default" and self.default_noise_std is None:
            raise ValueError("default_noise_std is required when unknown='default'")

        source_noise = {
            key: nonnegative_scalar(value, f"source_noise_std[{key!r}]")
            for key, value in self.source_noise_std.items()
        }
        default_noise = None
        if self.default_noise_std is not None:
            default_noise = nonnegative_scalar(self.default_noise_std, "default_noise_std")

        noise_values = []
        labels = []
        for source in source_values:
            if source in source_noise:
                noise = source_noise[source]
            elif normalized_unknown == "default":
                noise = default_noise
            else:
                raise ValueError(f"No noise_std was provided for source {source!r}")
            noise_values.append(noise)
            labels.append(f"{self.label_prefix}:{source}")

        return ObservationNoiseProfile(
            noise_std=np.asarray(noise_values, dtype=float),
            labels=np.asarray(labels, dtype=object),
            component_names=np.full(source_values.shape[0], "source", dtype=object),
        )

profile(sources)

Return per-row noise for source labels.

Source code in matgpr/noise_models.py
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
def profile(self, sources) -> ObservationNoiseProfile:
    """Return per-row noise for source labels."""
    source_values = _to_1d_object(sources, "sources")
    normalized_unknown = self.unknown.lower()
    if normalized_unknown not in {"error", "default"}:
        raise ValueError("unknown must be 'error' or 'default'")
    if normalized_unknown == "default" and self.default_noise_std is None:
        raise ValueError("default_noise_std is required when unknown='default'")

    source_noise = {
        key: nonnegative_scalar(value, f"source_noise_std[{key!r}]")
        for key, value in self.source_noise_std.items()
    }
    default_noise = None
    if self.default_noise_std is not None:
        default_noise = nonnegative_scalar(self.default_noise_std, "default_noise_std")

    noise_values = []
    labels = []
    for source in source_values:
        if source in source_noise:
            noise = source_noise[source]
        elif normalized_unknown == "default":
            noise = default_noise
        else:
            raise ValueError(f"No noise_std was provided for source {source!r}")
        noise_values.append(noise)
        labels.append(f"{self.label_prefix}:{source}")

    return ObservationNoiseProfile(
        noise_std=np.asarray(noise_values, dtype=float),
        labels=np.asarray(labels, dtype=object),
        component_names=np.full(source_values.shape[0], "source", dtype=object),
    )

ReplicateNoiseModel dataclass

Estimate observation noise from replicate target measurements.

Rows sharing the same group label are treated as replicate measurements. Groups with at least two observations receive their sample standard deviation. Singleton groups receive fallback_noise_std when supplied, otherwise the pooled replicate standard deviation, and finally min_noise_std if no replicate information is available.

Source code in matgpr/noise_models.py
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
210
211
212
213
214
215
216
217
218
219
220
221
222
223
@dataclass
class ReplicateNoiseModel:
    """Estimate observation noise from replicate target measurements.

    Rows sharing the same group label are treated as replicate measurements.
    Groups with at least two observations receive their sample standard
    deviation. Singleton groups receive ``fallback_noise_std`` when supplied,
    otherwise the pooled replicate standard deviation, and finally
    ``min_noise_std`` if no replicate information is available.
    """

    min_noise_std: float = 1e-8
    fallback_noise_std: float | None = None
    ddof: int = 1
    label_prefix: str = "replicate"

    def __post_init__(self) -> None:
        self.min_noise_std = nonnegative_scalar(self.min_noise_std, "min_noise_std")
        if self.fallback_noise_std is not None:
            self.fallback_noise_std = nonnegative_scalar(
                self.fallback_noise_std,
                "fallback_noise_std",
            )
        if self.ddof < 0:
            raise ValueError("ddof must be non-negative")

    def fit(self, y, groups):
        """Estimate group-level noise from targets and replicate labels."""
        y_values = to_1d_finite(y, "y")
        group_values = _to_1d_object(groups, "groups")
        validate_same_row_count(y_values, group_values, "y", "groups")

        frame = pd.DataFrame({"y": y_values, "group": group_values})
        group_noise: dict[object, float] = {}
        replicate_variances = []

        for group, group_frame in frame.groupby("group", sort=False, dropna=False):
            values = group_frame["y"].to_numpy(dtype=float)
            if values.shape[0] > self.ddof:
                noise = float(np.std(values, ddof=self.ddof))
                replicate_variances.append(noise**2)
            else:
                noise = np.nan
            group_noise[group] = noise

        if self.fallback_noise_std is not None:
            fallback = self.fallback_noise_std
        elif replicate_variances:
            fallback = float(np.sqrt(np.mean(replicate_variances)))
        else:
            fallback = self.min_noise_std

        self.group_noise_std_ = {
            group: max(noise if np.isfinite(noise) else fallback, self.min_noise_std)
            for group, noise in group_noise.items()
        }
        self.fallback_noise_std_ = max(fallback, self.min_noise_std)
        return self

    def profile(self, groups) -> ObservationNoiseProfile:
        """Return per-row noise for fitted replicate labels."""
        self._require_fitted()
        group_values = _to_1d_object(groups, "groups")
        noise_values = []
        labels = []
        for group in group_values:
            noise_values.append(self.group_noise_std_.get(group, self.fallback_noise_std_))
            labels.append(f"{self.label_prefix}:{group}")
        return ObservationNoiseProfile(
            noise_std=np.asarray(noise_values, dtype=float),
            labels=np.asarray(labels, dtype=object),
            component_names=np.full(group_values.shape[0], "replicate", dtype=object),
        )

    def fit_profile(self, y, groups) -> ObservationNoiseProfile:
        """Fit the model and return the training-row noise profile."""
        self.fit(y, groups)
        return self.profile(groups)

    def _require_fitted(self) -> None:
        if not hasattr(self, "group_noise_std_"):
            raise ValueError("ReplicateNoiseModel must be fitted before calling profile")

fit(y, groups)

Estimate group-level noise from targets and replicate labels.

Source code in matgpr/noise_models.py
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
def fit(self, y, groups):
    """Estimate group-level noise from targets and replicate labels."""
    y_values = to_1d_finite(y, "y")
    group_values = _to_1d_object(groups, "groups")
    validate_same_row_count(y_values, group_values, "y", "groups")

    frame = pd.DataFrame({"y": y_values, "group": group_values})
    group_noise: dict[object, float] = {}
    replicate_variances = []

    for group, group_frame in frame.groupby("group", sort=False, dropna=False):
        values = group_frame["y"].to_numpy(dtype=float)
        if values.shape[0] > self.ddof:
            noise = float(np.std(values, ddof=self.ddof))
            replicate_variances.append(noise**2)
        else:
            noise = np.nan
        group_noise[group] = noise

    if self.fallback_noise_std is not None:
        fallback = self.fallback_noise_std
    elif replicate_variances:
        fallback = float(np.sqrt(np.mean(replicate_variances)))
    else:
        fallback = self.min_noise_std

    self.group_noise_std_ = {
        group: max(noise if np.isfinite(noise) else fallback, self.min_noise_std)
        for group, noise in group_noise.items()
    }
    self.fallback_noise_std_ = max(fallback, self.min_noise_std)
    return self

profile(groups)

Return per-row noise for fitted replicate labels.

Source code in matgpr/noise_models.py
201
202
203
204
205
206
207
208
209
210
211
212
213
214
def profile(self, groups) -> ObservationNoiseProfile:
    """Return per-row noise for fitted replicate labels."""
    self._require_fitted()
    group_values = _to_1d_object(groups, "groups")
    noise_values = []
    labels = []
    for group in group_values:
        noise_values.append(self.group_noise_std_.get(group, self.fallback_noise_std_))
        labels.append(f"{self.label_prefix}:{group}")
    return ObservationNoiseProfile(
        noise_std=np.asarray(noise_values, dtype=float),
        labels=np.asarray(labels, dtype=object),
        component_names=np.full(group_values.shape[0], "replicate", dtype=object),
    )

fit_profile(y, groups)

Fit the model and return the training-row noise profile.

Source code in matgpr/noise_models.py
216
217
218
219
def fit_profile(self, y, groups) -> ObservationNoiseProfile:
    """Fit the model and return the training-row noise profile."""
    self.fit(y, groups)
    return self.profile(groups)

FeatureNoiseModel dataclass

Compute heteroscedastic noise from feature values.

The supplied callable receives a numeric feature matrix and returns one noise standard deviation per row. This lets users encode physics-aware assumptions, such as increasing measurement uncertainty at high temperature, low signal intensity, extreme compositions, or high load.

Source code in matgpr/noise_models.py
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
@dataclass(frozen=True)
class FeatureNoiseModel:
    """Compute heteroscedastic noise from feature values.

    The supplied callable receives a numeric feature matrix and returns one
    noise standard deviation per row. This lets users encode physics-aware
    assumptions, such as increasing measurement uncertainty at high
    temperature, low signal intensity, extreme compositions, or high load.
    """

    noise_std_function: Callable[[np.ndarray], Sequence[float] | np.ndarray]
    label: str = "feature_noise"

    def profile(self, X) -> ObservationNoiseProfile:
        """Return per-row noise from a feature-dependent equation."""
        x_values = to_2d_finite(_as_numeric_matrix(X), "X")
        noise = to_1d_finite(self.noise_std_function(x_values), "noise_std_function(X)")
        if noise.shape[0] != x_values.shape[0]:
            raise ValueError("noise_std_function must return one value per feature row")
        if np.any(noise < 0):
            raise ValueError("noise_std_function must return non-negative values")
        return ObservationNoiseProfile(
            noise_std=noise,
            labels=np.full(x_values.shape[0], self.label, dtype=object),
            component_names=np.full(x_values.shape[0], "feature", dtype=object),
        )

profile(X)

Return per-row noise from a feature-dependent equation.

Source code in matgpr/noise_models.py
239
240
241
242
243
244
245
246
247
248
249
250
251
def profile(self, X) -> ObservationNoiseProfile:
    """Return per-row noise from a feature-dependent equation."""
    x_values = to_2d_finite(_as_numeric_matrix(X), "X")
    noise = to_1d_finite(self.noise_std_function(x_values), "noise_std_function(X)")
    if noise.shape[0] != x_values.shape[0]:
        raise ValueError("noise_std_function must return one value per feature row")
    if np.any(noise < 0):
        raise ValueError("noise_std_function must return non-negative values")
    return ObservationNoiseProfile(
        noise_std=noise,
        labels=np.full(x_values.shape[0], self.label, dtype=object),
        component_names=np.full(x_values.shape[0], "feature", dtype=object),
    )

constant_noise_profile(n_observations, noise_std, *, label='constant')

Return a constant noise profile.

Source code in matgpr/noise_models.py
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
def constant_noise_profile(
    n_observations: int,
    noise_std: float,
    *,
    label: str = "constant",
) -> ObservationNoiseProfile:
    """Return a constant noise profile."""
    if n_observations <= 0:
        raise ValueError("n_observations must be positive")
    noise = np.full(int(n_observations), nonnegative_scalar(noise_std, "noise_std"), dtype=float)
    labels = np.full(int(n_observations), label, dtype=object)
    return ObservationNoiseProfile(
        noise_std=noise,
        labels=labels,
        component_names=np.full(int(n_observations), "constant", dtype=object),
    )

combine_noise_profiles(*profiles, mode='quadrature', label='combined')

Combine independent noise profiles.

mode="quadrature" combines independent variance components as sqrt(sum(sigma_i^2)). mode="max" keeps the largest standard deviation per row, and mode="sum" adds standard deviations directly.

Source code in matgpr/noise_models.py
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
def combine_noise_profiles(
    *profiles: ObservationNoiseProfile,
    mode: str = "quadrature",
    label: str = "combined",
) -> ObservationNoiseProfile:
    """Combine independent noise profiles.

    ``mode="quadrature"`` combines independent variance components as
    ``sqrt(sum(sigma_i^2))``. ``mode="max"`` keeps the largest standard
    deviation per row, and ``mode="sum"`` adds standard deviations directly.
    """
    if len(profiles) == 1 and isinstance(profiles[0], (list, tuple)):
        profiles = tuple(profiles[0])
    if not profiles:
        raise ValueError("At least one ObservationNoiseProfile is required")
    if not all(isinstance(profile, ObservationNoiseProfile) for profile in profiles):
        raise TypeError("profiles must contain only ObservationNoiseProfile objects")

    n_observations = profiles[0].n_observations
    for profile in profiles:
        if profile.n_observations != n_observations:
            raise ValueError("All noise profiles must have the same number of observations")

    noise_matrix = np.vstack([profile.noise_std for profile in profiles])
    normalized_mode = mode.lower()
    if normalized_mode == "quadrature":
        noise = np.sqrt(np.sum(noise_matrix * noise_matrix, axis=0))
    elif normalized_mode == "max":
        noise = np.max(noise_matrix, axis=0)
    elif normalized_mode == "sum":
        noise = np.sum(noise_matrix, axis=0)
    else:
        raise ValueError("mode must be one of: quadrature, max, sum")

    component_labels = np.asarray(
        [
            "+".join(str(profile.component_names[row]) for profile in profiles)
            for row in range(n_observations)
        ],
        dtype=object,
    )
    return ObservationNoiseProfile(
        noise_std=noise,
        labels=np.full(n_observations, label, dtype=object),
        component_names=component_labels,
    )