igmhub / LaCE

Emulator for the one-dimensional Lyman-alpha flux power spectrum. It was used to generate the results shown in https://arxiv.org/abs/2011.15127, https://arxiv.org/abs/2209.09895, and https://arxiv.org/abs/2305.19064.
1 stars 3 forks source link

Polyfit emulator covariance #3

Closed Chris-Pedersen closed 2 years ago

Chris-Pedersen commented 2 years ago

I noticed this was being calculated poorly when writing up the paper text. Here we try a better motivated treatment and find that it performs worse.

andreufont commented 2 years ago

I believe there are two problems with the derivation / code.

The first one is that GPy.predict returns an estimate of the variance, but we transform this to an error in the predict function in emulator_gp.py: out_err=np.ndarray.flatten(np.sqrt(err)*self.scalefactors)

However, in your derivation (and in the existing and new codes) you are treating err[n] as a variance. I believe you should square these terms.

The other issue is that the equation for the error should have all the terms log(k) squared, since the partial derivatives should be squared as well.

Taking into account both issues, I think the right code should be:

                    # first estimate error on y=log P
                    lk=np.log(k_Mpc)
                    erry2=((err[0]*lk**4)**2 + (err[1]*lk**3)**2
                            + (err[2]*lk**2)**2 + (err[3]*lk)**2 + err[4]**2)
                    # compute error on P
                    err=interpolated_P*np.sqrt(erry2)
andreufont commented 2 years ago

Using my suggested code above, I get more reasonable emulator errors. I attach two plots correpsonding to the ones in your notebook. new_test_1 new_test_2

andreufont commented 2 years ago

Having said that, the derivations assumes that the coefficients of the polynomials are independent, and this is probably not the case. It could be if we did a PCA decomposition or similar, but not with simple polynomials.

Anyway, I could re-run a couple of chains with this new estimate of the covariance (if you agree with my derivation), and then we can decide whether we use it or just ignore the emulator covariance alltogether.

Chris-Pedersen commented 2 years ago

Thanks for spotting both of these, had to do this in a bit of a rush last night. I think this fixes it and is the better procedure. In practice this looks like a similar emulator covariance to what we had before so the chains should come out looking the same, hope this doesn't require too much rerunning.

andreufont commented 2 years ago

Let's discuss these today with Nick. Babysitting the chains at NERSC might take some time, but that's ok, I'll do it if we have to.