Use fixed BayesEphem parameters by extending Parameter classes

[paulthebaker]

I've been thinking about using a fixed BayesEphem model for single pulsar noise runs. These are for model selection to determine what sort of RN model to use for that pulsar. I don't want an ephemeris error to be picked up so it seems like a more complicated model is needed. It may not matter, but I'd like to find out...

There are two options: 1) fix ephemeris params to mean values from full PTA GWB analysis 2) use informative priors on ephemeris params from full PTA GWB analysis output by fitting a normal distr to the chains

(this is only marginally cheating since most of the BE parameter info is coming from the other N-1 pulsars)

I can write a modified PhysicalEphemerisSignal class factory for each case... but I would prefer to extend the existing Parameter classes to handle the it. I may need some help.

Method 1

I want to do something like

eph = deterministic_signals.PhysicalEphemerisSignal(
                frame_drift_rate=parameter.Constant()('frame_drift_rate'),
                d_jupiter_mass=parameter.Constant()('d_jupiter_mass'),
                d_saturn_mass=parameter.Constant()('d_saturn_mass'),
                d_uranus_mass=parameter.Constant()('d_uranus_mass'),
                d_neptune_mass=parameter.Constant()('d_neptune_mass'),
                jup_orb_elements=parameter.Constant(size=6)('jup_orb_elements'),
                use_epoch_toas=True)

then use PTA.set_default_params() to populate them.

If we modify Constant to take a size input (like Normal and Uniform), then the above should just work.

Method 2

To do (2) we need to allow a size > 1 parameter to take a list of distribution parameters instead of one. Currently, we do things like this

jup_orb_elements=parameter.Uniform(-0.05, 0.05, size=6)('jup_orb_elements')

to get 6 orbital parameters with the same prior.

I would like to do this

joe_mean = [-0.0072, -0.0036, -0.0099, -0.0099, 0.0015, 0.0150]
joe_std = [0.0068, 0.0082, 0.0059, 0.0084, 0.0054, 0.0098]
jup_orb_elements=parameter.Normal(joe_mean, joe_std, size=6)('jup_orb_elements')

to get 6 parameters with different priors (although the same underlying distribution. The size could even be inferred from the size of the input arrays.

We could also allow

parameter.Normal([a, b, c], 1)

to get three parameters with different means, but the same standard deviation.

[paulthebaker]

This needs to wait until #119 is finished...

@vallis, I'm looking at the new Parameter in #119. It looks like

parameter.Normal([a, b, c], 1)

would create a 3D multivariate normal. Is this correct?

Instead we could pass a list of Parameters

joe_mean = [-0.0072, -0.0036, -0.0099, -0.0099, 0.0015, 0.0150]
joe_std = [0.0068, 0.0082, 0.0059, 0.0084, 0.0054, 0.0098]
jup_orb_elements=[parameter.Normal(mu, std)('jup_orb_elements') for mu, std in zip(joe_mean, joe_std)]

We would need to handle the multiple identically named parameters being properly broken out to jup_orb_elements_0, etc.