Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Normalised stats #396

Merged
merged 6 commits into from
Sep 27, 2024
Merged

Normalised stats #396

Show file tree
Hide file tree
Changes from 3 commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
7754bab
allow passing a `NestedSamples.stats` instance as a normalisation ref…
lukashergt Sep 19, 2024
fa31ee3
add tests for the normalisation kwarg `norm` for `NestedSamples.stats`
lukashergt Sep 19, 2024
44334bd
bump version to 2.9.0
lukashergt Sep 19, 2024
55a3928
improve docstring for `norm` kwarg, saying what columns will addition…
lukashergt Sep 26, 2024
60030df
Merge branch 'master' into normalised_stats
lukashergt Sep 27, 2024
fa919c1
Merge branch 'master' into normalised_stats
lukashergt Sep 27, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion README.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
anesthetic: nested sampling post-processing
===========================================
:Authors: Will Handley and Lukas Hergt
:Version: 2.8.14
:Version: 2.9.0
:Homepage: https://github.com/handley-lab/anesthetic
:Documentation: http://anesthetic.readthedocs.io/

Expand Down
2 changes: 1 addition & 1 deletion anesthetic/_version.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1 +1 @@
__version__ = '2.8.14'
__version__ = '2.9.0'
27 changes: 26 additions & 1 deletion anesthetic/samples.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -761,7 +761,7 @@ def ns_output(self, *args, **kwargs):
" as well as average loglikelihoods: help(samples.logL_P)"
)

def stats(self, nsamples=None, beta=None):
def stats(self, nsamples=None, beta=None, norm=None):
r"""Compute Nested Sampling statistics.

Using nested sampling we can compute:
Expand Down Expand Up @@ -822,6 +822,11 @@ def stats(self, nsamples=None, beta=None):
beta : float, array-like, optional
inverse temperature(s) beta=1/kT. Default self.beta

norm : Series, :class:`Samples`, optional
:meth:`NestedSamples.stats` output used for normalisation.
Can be either a Series of mean values or Samples produced with
matching `nsamples` and `beta`.

lukashergt marked this conversation as resolved.
Show resolved Hide resolved
Returns
-------
if beta is scalar and nsamples is None:
Expand Down Expand Up @@ -861,6 +866,26 @@ def stats(self, nsamples=None, beta=None):
samples.set_label('d_G', r'$d_\mathrm{G}$')

samples.label = self.label

if norm is not None:
lukashergt marked this conversation as resolved.
Show resolved Hide resolved
samples['Delta_logZ'] = samples['logZ'] - norm['logZ']
samples.set_label('Delta_logZ',
r"$\Delta\ln\mathcal{Z}$")

samples['Delta_D_KL'] = samples['D_KL'] - norm['D_KL']
samples.set_label('Delta_D_KL',
r"$\Delta\mathcal{D}_\mathrm{KL}$")

samples['Delta_logL_P'] = samples['logL_P'] - norm['logL_P']
samples.set_label(
'Delta_logL_P',
r"$\Delta\langle\ln\mathcal{L}\rangle_\mathcal{P}$"
)

samples['Delta_d_G'] = samples['d_G'] - norm['d_G']
samples.set_label('Delta_d_G',
r"$\Delta d_\mathrm{G}$")

return samples

def logX(self, nsamples=None):
Expand Down
54 changes: 54 additions & 0 deletions tests/test_samples.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -946,39 +946,93 @@ def test_stats():
beta = [0., 0.5, 1.]

vals = ['logZ', 'D_KL', 'logL_P', 'd_G']
delta_vals = ['Delta_logZ', 'Delta_D_KL', 'Delta_logL_P', 'Delta_d_G']

labels = [r'$\ln\mathcal{Z}$',
r'$\mathcal{D}_\mathrm{KL}$',
r'$\langle\ln\mathcal{L}\rangle_\mathcal{P}$',
r'$d_\mathrm{G}$']
delta_labels = [r'$\Delta\ln\mathcal{Z}$',
r'$\Delta\mathcal{D}_\mathrm{KL}$',
r'$\Delta\langle\ln\mathcal{L}\rangle_\mathcal{P}$',
r'$\Delta d_\mathrm{G}$']

stats = pc.stats()
assert isinstance(stats, WeightedLabelledSeries)
assert_array_equal(stats.drop_labels().index, vals)
assert_array_equal(stats.get_labels(), labels)

stats = pc.stats(norm=pc.stats())
assert isinstance(stats, WeightedLabelledSeries)
assert_array_equal(stats.drop_labels().index, vals + delta_vals)
assert_array_equal(stats.get_labels(), labels + delta_labels)

stats = pc.stats(nsamples=nsamples)
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals)
assert_array_equal(stats.get_labels(), labels)
assert stats.index.name == 'samples'
assert_array_equal(stats.index, range(nsamples))

stats = pc.stats(nsamples=nsamples, norm=pc.stats())
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals + delta_vals)
assert_array_equal(stats.get_labels(), labels + delta_labels)
assert stats.index.name == 'samples'
assert_array_equal(stats.index, range(nsamples))

stats = pc.stats(nsamples=nsamples, norm=pc.stats(nsamples=nsamples))
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals + delta_vals)
assert_array_equal(stats.get_labels(), labels + delta_labels)
assert stats.index.name == 'samples'
assert_array_equal(stats.index, range(nsamples))

stats = pc.stats(beta=beta)
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals)
assert_array_equal(stats.get_labels(), labels)
assert stats.index.name == 'beta'
assert_array_equal(stats.index, beta)

stats = pc.stats(beta=beta, norm=pc.stats())
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals + delta_vals)
assert_array_equal(stats.get_labels(), labels + delta_labels)
assert stats.index.name == 'beta'
assert_array_equal(stats.index, beta)

stats = pc.stats(beta=beta, norm=pc.stats(beta=beta))
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals + delta_vals)
assert_array_equal(stats.get_labels(), labels + delta_labels)
assert stats.index.name == 'beta'
assert_array_equal(stats.index, beta)

stats = pc.stats(nsamples=nsamples, beta=beta)
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals)
assert_array_equal(stats.get_labels(), labels)
assert stats.index.names == ['beta', 'samples']
assert stats.index.levshape == (len(beta), nsamples)

stats = pc.stats(nsamples=nsamples, beta=beta, norm=pc.stats())
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals + delta_vals)
assert_array_equal(stats.get_labels(), labels + delta_labels)
assert stats.index.names == ['beta', 'samples']
assert stats.index.levshape == (len(beta), nsamples)

stats = pc.stats(nsamples=nsamples, beta=beta,
norm=pc.stats(nsamples=nsamples, beta=beta))
assert isinstance(stats, WeightedLabelledDataFrame)
assert_array_equal(stats.drop_labels().columns, vals + delta_vals)
assert_array_equal(stats.get_labels(), labels + delta_labels)
assert stats.index.names == ['beta', 'samples']
assert stats.index.levshape == (len(beta), nsamples)

for beta in [1., 0., 0.5]:
np.random.seed(42)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why do you need a seed here?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I managed to get failing tests in what follows, because of changes to the seed from the tests I added above... Nothing dramatic, just bad chance in combination with the tests checking against 3*std or similar. I could have fiddled around with tolerances, but just trying a new seed there was enough, so I did not dig deeper...

Copy link
Collaborator

@AdamOrmondroyd AdamOrmondroyd Sep 27, 2024 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ok, I'll reapprove when you've done the version number ofc this is a minor bump

pc.beta = beta
n = 1000
PC = pc.stats(n, beta)
Expand Down
Loading