yabebalFantaye / MCEvidence

A python package implementing the MARGINAL LIKELIHOODS FROM MONTE CARLO MARKOV CHAINS algorithm described in Heavens et. al. (2017)
MIT License
38 stars 17 forks source link

Negative eigenvalues of the covariance matrix #12

Open subhajitghosh-phy opened 2 years ago

subhajitghosh-phy commented 2 years ago

Hello,

I am trying to use this code to compute the Bayesian Evidence for chains computed using MontePython. I tried to test the code first on official Planck chains.

However, I am getting the following error: " ValueError Traceback (most recent call last)

in () ----> 1 MCEvidence('/afs/crc.nd.edu/user/s/sghosh5/driso/planck-official-chains/base/plikHM_TT/base_plikHM_TT').evidence() /afs/crc.nd.edu/user/s/sghosh5/MCEvidence/MCEvidence.pyc in evidence(self, verbose, rand, info, covtype, profile, pvolume, pos_lnp, nproc, prewhiten) 1055 1056 #diagonalise samples to have unit variance -> 1057 samples = self.diagonalise_chain(samples, covstat['eVec'], covstat['eVal']) 1058 1059 #We need the logarithm of the likelihood - not the negative log /afs/crc.nd.edu/user/s/sghosh5/MCEvidence/MCEvidence.pyc in diagonalise_chain(self, s, eigenVec, eigenVal) 845 # And renormalise new parameters to have unit covariance matrix: 846 for i in range(self.ndim): --> 847 s[:,i]= s[:,i]/math.sqrt(eigenVal[i]) 848 849 return s ValueError: math domain error " I think the error is happening due to negative eigenvalues of the covariance matrix (the code is trying to evaluate the square root of the eigenvalue) -- since I got the following warning before that! " INFO:MCEvidence.py.__init__():781 chain array dimensions: [37489] x 84 = INFO:MCEvidence.py.get_covariance():859 Estimating covariance matrix using all chains INFO:MCEvidence.py.get_covariance():863 covariance matrix estimated using nsample=37489 WARNING:MCEvidence.py.get_covariance():869 Some of the eigenvalues of the covariance matrix are negative and/or complex: INFO:MCEvidence.py.get_samples():910 getting samples for partition s1: nsamples=37489 " Can anyone please tell me what to do here? The official chains I am using are from: http://pla.esac.esa.int/pla/aio/product-action?COSMOLOGY.FILE_ID=COM_CosmoParams_base-plikHM_R3.01.zip Thanks in advance! Best, Subhajit
subhajitghosh-phy commented 2 years ago

I think I figured it out! I looked into the associated Planck chain analysis code supplied with the code. As long as I restrict the dimensions to the primary cosmological parameters (by ndim), I am not getting negative eigenvalues and the evidence calculation goes through!

The negative eigenvalues all correspond to the nuisance parameters. I will be happy if someone can explain what negative eigenvalues means in this context!

Thanks!