feat: adding a samplers module with a basic ptmcmc wrapper script #50

stevertaylor · 2024-08-24T16:44:18Z

No description provided.

…r ptmcmc" This reverts commit ed16330.

AaronDJohnson

This could really save us time down the road with people asking us questions on how to use the models. I highly recommend that we do this sort of thing. Let's discuss the implementation a bit though.

AaronDJohnson · 2024-08-25T15:16:06Z

src/discovery/models.py

+                                                      common=['gw_log10_A']+gamma_common_name, name='gw')
+                                                      ]) for psr in psrs))
+    elif common_type == 'hd':
+        gl = ds.GlobalLikelihood((ds.PulsarLikelihood([psr.residuals,


See above. I think that this extra Likelihood block could be modified to fit within just one block if we make a small if then statement choosing between makegp_fourier and makegp_fourier_global.

I'm happy to do that. But I'm confused about ArrayLikelihood being adopted as the standard. The README mentions it as still being experimental. If we should adopt is as standard, then why does GlobalLikelihood exist at all?

AaronDJohnson · 2024-08-25T15:17:47Z

src/discovery/samplers/ptmcmc_wrapper.py

+from PTMCMCSampler import PTMCMCSampler as ptmcmc
+
+
+class JumpProposal(object):


These JumpProposals should really go in the sampler instead. @vhaasteren can we put them there instead?

I asked this below, but do you mean in PTMCMCSampler itself? We could do that, and it would certainly be a cleaner way for both discovery and enterprise to access these custom jump proposal schemes. The issue is that the sampler package then becomes super specific to PTA data analysis, and loses its generality. If we do go that route, then I'd advocate for the generality of PTMCMCSampler to be retained, but with some additional features that PTA data analysts can access without being intrusive or obvious to the more general user.

AaronDJohnson · 2024-08-25T15:24:09Z

src/discovery/models.py

+from jax.tree_util import Partial
+import discovery as ds
+
+def lhood_maker(psrs, noisedict=None, gamma_common=None, 


I think make_likelihood would be a better name for this.

AaronDJohnson · 2024-08-25T15:24:24Z

src/discovery/models.py

+        gamma_common_name = ['gw_gamma']
+
+    if common_type == 'curn':
+        gl = ds.GlobalLikelihood((ds.PulsarLikelihood([psr.residuals,


We should exclusively use ArrayLikelihood in these functions. It is faster, and these models typically assume that the noise is the "same." If they do not, these e_e type functions won't save us and the full GlobalLikelihood functionality needs to be visited.

AaronDJohnson · 2024-08-25T15:26:21Z

src/discovery/samplers/ptmcmc_wrapper.py

+
+class JumpProposal(object):
+
+    def __init__(self, param_names, empirical_distr=None,


empirical_distr isn't used here

AaronDJohnson · 2024-08-25T15:27:06Z

src/discovery/samplers/ptmcmc_wrapper.py

+class JumpProposal(object):
+
+    def __init__(self, param_names, empirical_distr=None,
+                 save_ext_dists=False, outdir='./chains'):


save_ext_dists and outdir aren't used here

Not yet. It's basically ported over from e_e, anticipating when we may want to use empirical distributions as proposals.

AaronDJohnson · 2024-08-25T15:34:17Z

src/discovery/samplers/ptmcmc_wrapper.py

+        self.get_lnlikelihood = lambda x: float(jlogl(self.x2p(x)))
+        self.get_lnprior = lambda x: float(jlogp(self.x2p(x)))
+
+    def x2p(self, x):


x2p should be defined outside of this class, because all samplers may have use for it.

AaronDJohnson · 2024-08-25T15:34:36Z

src/discovery/samplers/ptmcmc_wrapper.py

+        # does not handle vector parameters
+        return {par: val for par, val in zip(self.param_names, x)}
+
+    def p2x(self, p):


p2x should be defined outside of this class, because all samplers may have use for it.

AaronDJohnson · 2024-08-25T15:47:14Z

src/discovery/samplers/ptmcmc_wrapper.py

+        return q, float(lqxy)
+
+
+class InferenceModel(object):


I'm not sold on this class being necessary. It seems to me that this could just be a function (setup_ptmcmc_sampler or similar) that returns the sampler object.

There's more to this class though. The initialization takes in the discovery likelihood and the pulsar names, then prepares things for PTMCMC (like get_lnlikelihood and get_lnprior) in a way that feels familiar to enterprise users. Plus you have the separate function that prepares parameter groupings. We may be able to split some of this up, but we do need most of it.

AaronDJohnson · 2024-08-25T15:51:06Z

src/discovery/samplers/ptmcmc_wrapper.py

+        """
+        return np.array(list(p.values()), 'd')
+
+    def get_parameter_groups(self):


This function also seems like something that should be in the sampler package since it gets used in multiple places (here and e_e).

Do you mean in the PTMCMCSampler package itself? If so, that's certainly possible, but it would then make that package highly specific to PTA data analysis in a way that it is currently not, and would lose some of its generality as a sampling package. We could go that route though, but it's a bigger conversation.

stevertaylor added 4 commits August 23, 2024 20:27

adding new models module

5d2967a

feat: adding a samplers module with a basic wrapper script for ptmcmc

ed16330

Revert "feat: adding a samplers module with a basic wrapper script fo…

7aed3c0

…r ptmcmc" This reverts commit ed16330.

feat: adding new samplers module with basic ptmcmc wrapper script

1d1efb5

AaronDJohnson requested changes Aug 25, 2024

View reviewed changes

moving p2x and x2p into utils script inside samplers module

a183bf7

stevertaylor requested a review from AaronDJohnson September 16, 2024 17:09

		from PTMCMCSampler import PTMCMCSampler as ptmcmc


		class JumpProposal(object):


		class JumpProposal(object):

		def __init__(self, param_names, empirical_distr=None,

feat: adding a samplers module with a basic ptmcmc wrapper script #50

Are you sure you want to change the base?

feat: adding a samplers module with a basic ptmcmc wrapper script #50

Uh oh!

Conversation

stevertaylor commented Aug 24, 2024

Uh oh!

AaronDJohnson left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants