BellaNasirudin / py21cmmc_fg

A foreground and instrumental plugin for py21cmmc
MIT License
1 stars 4 forks source link

Gridded Centres, Thermal Noise, No Foregrounds #21

Open steven-murray opened 5 years ago

steven-murray commented 5 years ago

Here is the plots I get from running gridded_noise_nofg.py with DEBUG=False, under commit 7f543c0939c61a1d5dc5f3fb773999baeefdf2c0.

The plots.py script should be updated to get nicer plots, but i just wanted to get a rough picture.

Parameters of relevance (look at commit, especially base_definitions.py, for details): HII_DIM=250, BOX_LEN=750 Mpc. Max bl length=300.0., blackman taper.

instrumentalgridtestnoise_traceplot

Mock Power Spectrum: instrumentalgridtestnoise_dataps

Signal-to-noise: instrumentalgridtestnoise_signalnoiseps

Covariance at lowest u: instrumentalgridtestnoise_noisecov

Clearly not converging, in fact, each walker is having a lot of trouble finding any other position which has better likelihood. Acceptance fraction is lower than 10%. Will investigate more.

steven-murray commented 5 years ago

For commit a0e60bf94ffcc77049b430c4743d1d5f7e834e9e --

instrumentalgridtestnoise_traceplot instrumentalgridtestnoise_dataps instrumentalgridtestnoise_noiseps instrumentalgridtestnoise_noisecov

Clearly, things are going very wrong now.

Positively, the noise properties look reasonable... it is purely a function of |u|, and decreases linearly with |u| (since we average in circles). The covariance is also diagonal, which is true when there's no foregrounds. The thermal noise looks very low, if anything. So really I am unsure why the likelihood is acting so strangely.

steven-murray commented 5 years ago

Updated results from 949ec9d0cc8d4f4b479181a2cb0290753ee8af84. Run with "python gridded_noise_nofg.py", i.e. with DEBUG set to False (or zero).

instrumentalgridtestnoise_traceplot instrumentalgridtestnoise_2dps instrumentalgridtestnoise_weighting

It is still giving the huge bias (cf. first plot, the traceplot), which is really starting to frustrate me.

The second plot here shows a compendium of information from the output 2D PS. It shows the measured "mock" 2D power spectrum (i.e. after all baselines and thermal noise have been added in). Then it shows two "model" 2D PS, the first in the chain and the last. Recall that the model PS have the mean noise added, but not a random noise variable. The mean should look close to the model (within the std. dev.), which it does at a glance. One possible oddity is that the first and last 2D PS are very similar, even though their parameters are quite different. This could be a bug in the way that the parameters are (not) updated in the code, or it could be a high level of covariance between the parameters, which may explain the bias in the presence of noise. The lower panels of the second plot show noise properties. The left-most shows the model standard deviation expected in each cell, from thermal noise alone. It seems to have the correct form, which is to say it decreases linearly with increasing u. The right-most plots show how many sigma away from the data each model is (i.e. the model above that particular panel). The negative of these numbers are added to get the actual likelihood. Again, the only real oddity that jumps out is that they are so similar.

Finally, the third plot shows the baseline weights. The left shows the number of baselines in each cell (in the first frequency bin). The middle shows the effective number of baselines (wieghted by the taper and channel width, squared, cf the notebook on this in the devel/ directory) in each UV cell. The right show this effective number, but as a function of |u|, i.e. the final 2D PS perpendicular bins. It also shows the number of UV cells that go into each |u| bin. This all looks reasonable to me, so it doesn't seem like anything is going awry at this end.

I have only two ideas at this point: (i) somehow the parameters of the EoR are not being updated properly, so that very different parameters give similar or the same results. (ii) the analytic thermal noise model is incorrect somehow (either it is derived incorrectly, or not being numerically applied correctly).

The first can be tracked by calculating as extra data the 2D PS of the raw lightcone, using no baselines. In fact, we could even store slices of the lightcone itself to convince ourselves that the parameters are doing what we think they should.

The second can be formally tested by literally running a bunch of thermal noise simulations and checking the output mean/var against the analytic computation. I will try this first, and also run a test with the latter.

I am also going to increase the tile diameter to 20m so that not so much tiling needs to be done, so that smaller simulations can be used for this testing, without loss of physical resolution.

steven-murray commented 5 years ago

Starting to get better...

So, having had a look at point (ii) above (i.e. checking if analytic thermal noise model is incorrect somehow) showed that indeed it was. There were two factors of two that needed to be changed in different places. They have now been fixed.

I now get the following plot:

thermal_noise_test

This shows the thermal noise power mean and variance, for a numerical realisation, and analytic derivation. They are roughly the same, though there is definitely a consistent offset of the mean by 95% or so. I am not sure where this comes from.

Having re-run the MCMC chain, I get the following:

instrumentalgridtestnoise_2dps instrumentalgridtestnoise_cornerplot instrumentalgridtestnoise_traceplot

This was supposed to just be a short run to see if the chains wouldn't run off somewhere they shouldn't, and as such it hasn't fully converged. But we can see already that in fact it is running off somewhere it shouldn't. At least now the bias is in the other direction (actually the same direction the bias goes in if we use no noise). The bias is still very large, but arguably smaller than before. This might be because of that 5% discrepancy, but I am unsure (and if it is, not sure where that comes from or how to fix it).

I also added tracking of a slice of the raw lightcone, and I found that the parameters definitely do update the field, though most of that seems to be in the mean brightness, which may be one reason that it is not converging well.

I also just ran a version without thermal noise (I haven't done this in a while, and lots has been updated in the meantime, so I'm checking that it still works, roughly) which seems to work fine:

instrumentalgridtest_traceplot

and I'm also running a version where the thermal noise is calculated numerically, which should help check whether the problem is in that 5% discrepancy.

Other than that, I am still very unsure. We could try expanding the frequency range, and/or making the tile diameter much smaller (i.e. not tiling as much, so we can access smaller perpendicular scales).

steven-murray commented 5 years ago

So, it would seem that the remaining problem is the 5% discrepancy. Running with numerically-calculated estimates of thermal noise power gives the following:

instrumentalgridtestnoise_2dps instrumentalgridtestnoise_cornerplot instrumentalgridtestnoise_traceplot

The traceplot seems to be converging reasonably. The lower right panels of the 2D PS plots indicate that the mock is more symmetric around the model prediction, which is what we should expect. So we have two options: (i) continue by using numerically-produced estimates of the noise (both foregrounds and thermal) (ii) try to figure out why the thermal noise is wrong.

At the moment there is no simple mechanism to force the likelihood to use the numerical estimate of thermal noise (though you can just add in a foreground core that adds zero foregrounds to the sky to trigger it). This is probably worth adding in at some point, but at this stage it is probably better to try to forge ahead as quickly as possible.

Thus I suggest moving to tests that include point-source foregrounds, using the simple numerical estimate of mean/cov. We can come back and tidy up as we see fit.

I am now running a bigger version of this chain, to assess it in more detail.

steven-murray commented 5 years ago

Update: running longer chain with bigger boxes (i.e. DEBUG=0, commit 6faf91d0c8c62d5f22b7820c62e0c0b3e48996a6) with NUMERICAL estimate of thermal noise gives:

instrumentalgridtestnoisenumerical_2dps instrumentalgridtestnoisenumerical_cornerplot instrumentalgridtestnoisenumerical_traceplot

It looks like everything is working just fine. I am trying to work through why the analytic prediction does not match the numerical one. Will report back ASAP.