Skip to content

Commit 937137c

Browse files
committed
feat(families): score for builtins is implemented
1 parent 29019d1 commit 937137c

9 files changed

Lines changed: 611 additions & 2 deletions

File tree

src/pysatl_core/families/builtins/continuous/exponential.py

Lines changed: 56 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -223,6 +223,38 @@ def _support(_: Parametrization) -> ContinuousSupport:
223223
"""Support of exponential distribution"""
224224
return ContinuousSupport(left=0.0)
225225

226+
def _base_score(parameters: Parametrization, x: NumericArray) -> NumericArray:
227+
"""
228+
Compute the score (gradient of log‑PDF) for the base parametrization (rate λ).
229+
230+
For exponential distribution with rate λ > 0, log‑PDF is:
231+
log f(x) = log λ - λ x for x ≥ 0, else -∞.
232+
233+
The derivative with respect to λ is:
234+
∂/∂λ log f = 1/λ - x (for x ≥ 0).
235+
236+
For points x < 0 the density is zero; we return 0 for numerical stability
237+
(though the score is technically undefined there).
238+
239+
Parameters
240+
----------
241+
parameters : Parametrization
242+
Base parametrization instance (_Rate) with field lambda_.
243+
x : NumericArray
244+
Points at which to evaluate the gradient.
245+
246+
Returns
247+
-------
248+
NumericArray
249+
Gradient array of shape (..., 1) where the last axis corresponds to
250+
d(log f)/d(lambda_).
251+
"""
252+
params = cast(_Rate, parameters)
253+
lam = params.lambda_
254+
inside = x >= 0
255+
grad = np.where(inside, 1.0 / lam - x, 0.0)
256+
return grad[..., np.newaxis] # shape (..., 1)
257+
226258
Exponential = ParametricFamily(
227259
name=FamilyName.EXPONENTIAL,
228260
distr_type=UnivariateContinuous,
@@ -239,6 +271,7 @@ def _support(_: Parametrization) -> ContinuousSupport:
239271
CharacteristicName.KURT: kurt_func,
240272
},
241273
support_by_parametrization=_support,
274+
base_score=_base_score,
242275
)
243276
Exponential.__doc__ = EXPONENTIAL_DOC
244277

@@ -289,4 +322,27 @@ def transform_to_base_parametrization(self) -> Parametrization:
289322
"""
290323
return _Rate(lambda_=1.0 / self.beta)
291324

325+
def gradient_transform(self, grad_base: NumericArray) -> NumericArray:
326+
"""
327+
Transform gradient from base parameter (rate λ) to scale β = 1/λ.
328+
329+
The transformation is: β = 1/λ.
330+
Derivative: dβ/dλ = -1/λ².
331+
Hence: grad_β = grad_λ * (-1/λ²) = -grad_λ / λ².
332+
333+
Parameters
334+
----------
335+
grad_base : NumericArray
336+
Gradient with respect to λ, shape (..., 1).
337+
338+
Returns
339+
-------
340+
NumericArray
341+
Gradient with respect to β, shape (..., 1).
342+
"""
343+
lam = 1.0 / self.beta # λ = 1/β
344+
grad_lam = grad_base[..., 0]
345+
grad_beta = -grad_lam / (lam * lam)
346+
return grad_beta[..., np.newaxis].astype(np.float64)
347+
292348
ParametricFamilyRegister.register(Exponential)

src/pysatl_core/families/builtins/continuous/normal.py

Lines changed: 110 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -232,6 +232,36 @@ def _support(_: Parametrization) -> ContinuousSupport:
232232
"""Support of normal distribution"""
233233
return ContinuousSupport()
234234

235+
def _base_score(parameters: Parametrization, x: NumericArray) -> NumericArray:
236+
"""
237+
Compute the score (gradient of log‑PDF) for the base parametrization.
238+
239+
The base parametrization uses parameters mu (mean) and sigma (standard deviation).
240+
The returned gradient has shape (..., 2) with the last axis corresponding
241+
to [d(log f)/d(mu), d(log f)/d(sigma)].
242+
243+
Parameters
244+
----------
245+
parameters : Parametrization
246+
Base parametrization instance (must be _MeanStd).
247+
x : NumericArray
248+
Points at which the gradient is evaluated.
249+
250+
Returns
251+
-------
252+
NumericArray
253+
Gradient array of shape (..., 2).
254+
"""
255+
params = cast(_MeanStd, parameters)
256+
mu = params.mu
257+
sigma = params.sigma
258+
259+
z = (x - mu) / sigma
260+
261+
grad_mu = z / sigma
262+
grad_sigma = (z * z - 1) / sigma
263+
return np.stack([grad_mu, grad_sigma], axis=-1)
264+
235265
Normal = ParametricFamily(
236266
name=FamilyName.NORMAL,
237267
distr_type=UnivariateContinuous,
@@ -248,6 +278,7 @@ def _support(_: Parametrization) -> ContinuousSupport:
248278
CharacteristicName.KURT: kurt_func,
249279
},
250280
support_by_parametrization=_support,
281+
base_score=_base_score,
251282
)
252283
Normal.__doc__ = NORMAL_DOC
253284

@@ -305,6 +336,29 @@ def transform_to_base_parametrization(self) -> Parametrization:
305336
sigma = math.sqrt(self.var)
306337
return _MeanStd(mu=self.mu, sigma=sigma)
307338

339+
def gradient_transform(self, grad_base: NumericArray) -> NumericArray:
340+
"""
341+
Transform gradient from base parameters (mu, sigma) to (mu, var).
342+
343+
The transformation is: var = sigma^2 → d(var) = 2 sigma d(sigma).
344+
345+
Parameters
346+
----------
347+
grad_base : NumericArray
348+
Gradient with respect to (mu, sigma), shape (..., 2).
349+
350+
Returns
351+
-------
352+
NumericArray
353+
Gradient with respect to (mu, var), shape (..., 2).
354+
"""
355+
sigma = np.sqrt(self.var)
356+
357+
grad_mu = grad_base[..., 0]
358+
grad_sigma = grad_base[..., 1]
359+
grad_var = grad_sigma / (2 * sigma)
360+
return np.stack([grad_mu, grad_var], axis=-1)
361+
308362
@parametrization(family=Normal, name="meanPrec")
309363
class _MeanPrec(Parametrization):
310364
"""
@@ -338,6 +392,30 @@ def transform_to_base_parametrization(self) -> Parametrization:
338392
sigma = math.sqrt(1 / self.tau)
339393
return _MeanStd(mu=self.mu, sigma=sigma)
340394

395+
def gradient_transform(self, grad_base: NumericArray) -> NumericArray:
396+
"""
397+
Transform gradient from base parameters (mu, sigma) to (mu, tau).
398+
399+
The transformation is: tau = 1/sigma^2 → d(tau) = -2 / sigma^3 d(sigma).
400+
401+
Parameters
402+
----------
403+
grad_base : NumericArray
404+
Gradient with respect to (mu, sigma), shape (..., 2).
405+
406+
Returns
407+
-------
408+
NumericArray
409+
Gradient with respect to (mu, tau), shape (..., 2).
410+
"""
411+
412+
sigma = 1.0 / np.sqrt(self.tau)
413+
414+
grad_mu = grad_base[..., 0]
415+
grad_sigma = grad_base[..., 1]
416+
grad_tau = -grad_sigma * sigma**3 / 2.0
417+
return np.stack([grad_mu, grad_tau], axis=-1)
418+
341419
@parametrization(family=Normal, name="exponential")
342420
class _Exp(Parametrization):
343421
"""
@@ -384,4 +462,36 @@ def transform_to_base_parametrization(self) -> Parametrization:
384462
sigma = math.sqrt(-1 / (2 * self.a))
385463
return _MeanStd(mu=mu, sigma=sigma)
386464

465+
def gradient_transform(self, grad_base: NumericArray) -> NumericArray:
466+
"""
467+
Transform gradient from base parameters (mu, sigma) to (a, b).
468+
469+
The transformation is:
470+
mu = -b/(2a), sigma = sqrt(-1/(2a)).
471+
The Jacobian is applied via the chain rule.
472+
473+
Parameters
474+
----------
475+
grad_base : NumericArray
476+
Gradient with respect to (mu, sigma), shape (..., 2).
477+
478+
Returns
479+
-------
480+
NumericArray
481+
Gradient with respect to (a, b), shape (..., 2).
482+
"""
483+
a = self.a
484+
b = self.b
485+
486+
dmu_da = b / (2 * a * a)
487+
dmu_db = -1 / (2 * a)
488+
dsigma_da = 1 / (2 * np.sqrt(-2 * a**3))
489+
dsigma_db = 0.0
490+
491+
grad_mu = grad_base[..., 0]
492+
grad_sigma = grad_base[..., 1]
493+
grad_a = grad_mu * dmu_da + grad_sigma * dsigma_da
494+
grad_b = grad_mu * dmu_db + grad_sigma * dsigma_db
495+
return np.stack([grad_a, grad_b], axis=-1)
496+
387497
ParametricFamilyRegister.register(Normal)

src/pysatl_core/families/builtins/continuous/uniform.py

Lines changed: 103 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -260,6 +260,44 @@ def _support(parameters: Parametrization) -> ContinuousSupport:
260260
right_closed=True,
261261
)
262262

263+
def _base_score(parameters: Parametrization, x: NumericArray) -> NumericArray:
264+
"""
265+
Compute the score (gradient of log‑PDF) for the base
266+
parametrization (lower_bound, upper_bound).
267+
268+
For uniform distribution on [a, b] with a < b, log‑PDF is:
269+
log f(x) = -log(b - a) for a ≤ x ≤ b, else -∞.
270+
271+
The gradient with respect to a and b (where defined, i.e., inside the support) is:
272+
∂/∂a log f = 1/(b - a)
273+
∂/∂b log f = -1/(b - a)
274+
275+
For points outside the support, the gradient is set to 0 (since density is zero,
276+
but the score is typically considered undefined; we return 0 for numerical safety).
277+
278+
Parameters
279+
----------
280+
parameters : Parametrization
281+
Base parametrization instance (_Standard) with fields lower_bound, upper_bound.
282+
x : NumericArray
283+
Points at which to evaluate the gradient.
284+
285+
Returns
286+
-------
287+
NumericArray
288+
Gradient array of shape (..., 2) where last axis corresponds to
289+
[d(log f)/d(lower_bound), d(log f)/d(upper_bound)].
290+
"""
291+
params = cast(_Standard, parameters)
292+
a = params.lower_bound
293+
b = params.upper_bound
294+
width = b - a
295+
296+
inside = (x >= a) & (x <= b)
297+
grad_a = np.where(inside, 1.0 / width, 0.0)
298+
grad_b = np.where(inside, -1.0 / width, 0.0)
299+
return np.stack([grad_a, grad_b], axis=-1)
300+
263301
Uniform = ParametricFamily(
264302
name=FamilyName.CONTINUOUS_UNIFORM,
265303
distr_type=UnivariateContinuous,
@@ -276,6 +314,7 @@ def _support(parameters: Parametrization) -> ContinuousSupport:
276314
CharacteristicName.KURT: kurt_func,
277315
},
278316
support_by_parametrization=_support,
317+
base_score=_base_score,
279318
)
280319
Uniform.__doc__ = UNIFORM_DOC
281320

@@ -333,6 +372,38 @@ def transform_to_base_parametrization(self) -> Parametrization:
333372
half_width = self.width / 2
334373
return _Standard(lower_bound=self.mean - half_width, upper_bound=self.mean + half_width)
335374

375+
def gradient_transform(self, grad_base: NumericArray) -> NumericArray:
376+
"""
377+
Transform gradient from base parameters (a, b) to (mean, width).
378+
379+
The transformation is:
380+
mean = (a + b) / 2
381+
width = b - a
382+
383+
The Jacobian matrix:
384+
d(mean)/da = 0.5, d(mean)/db = 0.5
385+
d(width)/da = -1, d(width)/db = 1
386+
387+
Hence:
388+
grad_mean = 0.5 * (grad_a + grad_b)
389+
grad_width = -grad_a + grad_b
390+
391+
Parameters
392+
----------
393+
grad_base : NumericArray
394+
Gradient with respect to (a, b), shape (..., 2).
395+
396+
Returns
397+
-------
398+
NumericArray
399+
Gradient with respect to (mean, width), shape (..., 2).
400+
"""
401+
grad_a = grad_base[..., 0]
402+
grad_b = grad_base[..., 1]
403+
grad_mean = 0.5 * (grad_a + grad_b)
404+
grad_width = -grad_a + grad_b
405+
return np.stack([grad_mean, grad_width], axis=-1)
406+
336407
@parametrization(family=Uniform, name="minRange")
337408
class _MinRange(Parametrization):
338409
"""
@@ -365,4 +436,36 @@ def transform_to_base_parametrization(self) -> Parametrization:
365436
"""
366437
return _Standard(lower_bound=self.minimum, upper_bound=self.minimum + self.range_val)
367438

439+
def gradient_transform(self, grad_base: NumericArray) -> NumericArray:
440+
"""
441+
Transform gradient from base parameters (a, b) to (minimum, range_val).
442+
443+
The transformation is:
444+
minimum = a
445+
range_val = b - a
446+
447+
The Jacobian matrix:
448+
d(minimum)/da = 1, d(minimum)/db = 0
449+
d(range_val)/da = -1, d(range_val)/db = 1
450+
451+
Hence:
452+
grad_minimum = grad_a
453+
grad_range = -grad_a + grad_b
454+
455+
Parameters
456+
----------
457+
grad_base : NumericArray
458+
Gradient with respect to (a, b), shape (..., 2).
459+
460+
Returns
461+
-------
462+
NumericArray
463+
Gradient with respect to (minimum, range_val), shape (..., 2).
464+
"""
465+
grad_a = grad_base[..., 0]
466+
grad_b = grad_base[..., 1]
467+
grad_minimum = grad_a
468+
grad_range = -grad_a + grad_b
469+
return np.stack([grad_minimum, grad_range], axis=-1)
470+
368471
ParametricFamilyRegister.register(Uniform)

0 commit comments

Comments
 (0)