search

brinckmann / montepython_public

Public repository for the Monte Python Code
MIT License
93 stars 77 forks source link

R-1=0 : Zero Gelman-Rubin convergence criteria #311

Closed Amlan1996 closed 9 months ago

Amlan1996 commented 1 year ago

Hi,

I am trying to run montepython to get the best-fit parameters for a new model with two extra parameters from 'scf_parameters' in addition to the six Lambda CDM parameters. Even though my class (modified for the above model) runs and gives the output well within the Lambda CDM range for a range of different values of the two scf_parametrs, the Gelman-Rubin convergence criteria R-1 gives zero for all the parameters. As a result, the covariance matrix is not being computed, and I cannot get any decent output of the posterior.

I also used the values given in the first column of Table 2 in the Planck paper: arXiv:1807.06209 to check whether my class is showing the correct output and for the input of 'H0 = 66.88', 'z_reio = 7.5', 'omega_b = 0.02212', 'omega_dm = 0.1206' (acting as CDM in my model), 'Omega_fld = 0',' Omega_scf = -1',' Omega_Lambda = 0', and 'scf_parameters = 0.0, 1e-4, 1.0, 0.0, 1e-4, 0.1, 0.0', I got the outputs within the valid range given paper except for recombination redshift which is coming around 1089.10 which is a bit outside the range, but I don't believe just for that I am getting the zero R-1 value for all the parameters. I tried all possible values of the two scf_paramters, attempted to increase their width, and tried to reduce sigma, but nothing worked.

Here is the glimpse of the .param file :

data.experiments=['Planck_highl_TTTEEE', 'Planck_lowl_EE', 'Planck_lowl_TT', 'S8_KIDS1000']

data.over_sampling=[1, 5]

# Cosmological parameters list

data.parameters['omega_b']      = [  2.2377,   None, None,      0.015, 0.01, 'cosmo']
data.parameters['omega_dm']    = [ 0.12010,   None, None,     0.0013,    1, 'cosmo']
data.parameters['H0']           = [ 66.88,   None, None,    0.1,    1, 'cosmo']
data.parameters['ln10^{10}A_s'] = [  3.0447,   None, None,      0.015,    1, 'cosmo']
data.parameters['n_s']          = [  0.9659,   None, None,     0.0042,    1, 'cosmo']
data.parameters['z_reio']      = [ 7.0,  3.0, 15.0,  1.0,    1, 'cosmo']

data.parameters['scf_parameters__1'] = [0.0,    None,  None,  0.0, 1, 'cosmo']
data.parameters['scf_parameters__2'] = [0.000055,  1e-5,  1e-4, 1e-6, 1, 'cosmo']
data.parameters['scf_parameters__3'] = [1.0,    None, None, 0.0, 1, 'cosmo']
data.parameters['scf_parameters__4'] = [0.0,    None, None, 0.0, 1, 'cosmo']
data.parameters['scf_parameters__5'] = [0.000055,  1e-5, 1e-4, 1e-6, 1, 'cosmo']
data.parameters['scf_parameters__6'] = [0.1,   None, None, 0.0, 1, 'cosmo']
data.parameters['scf_parameters__7'] = [0.0,    None, None, 0.0, 1, 'cosmo']

The nuisance parameters are the same as in the input/base2018TTTEEE.param

# Derived parameters

data.parameters['tau_reio']          = [1, None, None, 0,     1,   'derived']
data.parameters['Omega0_scf']       = [1, None, None, 0,     1,   'derived']
data.parameters['YHe']             = [1, None, None, 0,     1,   'derived']
#data.parameters['H0']              = [0, None, None, 0,     1,   'derived']
data.parameters['A_s']             = [0, None, None, 0,  1e-9,   'derived']
data.parameters['sigma8']          = [0, None, None, 0,     1,   'derived']
data.parameters['100*theta_s']          = [1, None, None, 0,     1,   'derived']

# Other cosmo parameters (fixed parameters, precision parameters, etc.)

data.cosmo_arguments['sBBN file'] = data.path['cosmo']+'/external/bbn/sBBN.dat'
# BBN file path is automatically set to match CLASS version if 'sBBN file' is requested
# You can force the code to use the exact BBN file passed above with flag
#data.custom_bbn_file = True

data.cosmo_arguments['k_pivot'] = 0.05

# The base model features two massless
# and one massive neutrino with m=0.06eV.
# The settings below ensures that Neff=3.046
# and m/omega = 93.14 eV
data.cosmo_arguments['N_ur'] = 3.044
data.cosmo_arguments['N_ncdm'] = 1
data.cosmo_arguments['m_ncdm'] = 0.06
data.cosmo_arguments['T_ncdm'] = 0.71611
data.cosmo_arguments['omega_cdm'] = 1.e-9
data.cosmo_arguments['Omega_Lambda'] = 0
data.cosmo_arguments['Omega_dcdmdr'] = 0
data.cosmo_arguments['Omega_k'] = 0
data.cosmo_arguments['T_cmb'] = 2.7255
data.cosmo_arguments['Omega_idm'] = 0
data.cosmo_arguments['Omega_fld'] = 0
data.cosmo_arguments['Omega_scf'] = -1e-5
data.cosmo_arguments['attractor_ic_scf'] = 'no'

# These two are required to get sigma8 as a derived parameter
# (class must compute the P(k) until sufficient k)
data.cosmo_arguments['output'] = 'mPk'
data.cosmo_arguments['P_k_max_h/Mpc'] = 1.

#------ Mcmc parameters ----

data.N=10
data.write_step=5

I am getting no errors while running this file and yet in the output file I am getting the R-1 values like this 

Computing mean values
--> Computing variance
--> Computing convergence criterium (Gelman-Rubin)
 -> R-1 is 0.000000     for  omega_b
           0.000000     for  omega_dm
           0.000000     for  H0
           0.000000     for  ln10^{10}A_s
           0.000000     for  n_s
           0.000000     for  z_reio
           0.000000     for  scf_parameters__2
           0.000000     for  scf_parameters__5
           0.000000     for  A_cib_217
           0.000000     for  xi_sz_cib
           0.000000     for  A_sz
           0.000000     for  ps_A_100_100
           0.000000     for  ps_A_143_143
           0.000000     for  ps_A_143_217
           0.000000     for  ps_A_217_217
           0.000000     for  ksz_norm
           0.000000     for  gal545_A_100
           0.000000     for  gal545_A_143
           0.000000     for  gal545_A_143_217
           0.000000     for  gal545_A_217
           0.000000     for  galf_TE_A_100
           0.000000     for  galf_TE_A_100_143
           0.000000     for  galf_TE_A_100_217
           0.000000     for  galf_TE_A_143
           0.000000     for  galf_TE_A_143_217
           0.000000     for  galf_TE_A_217
           0.000000     for  calib_100T
           0.000000     for  calib_217T
           0.000000     for  A_planck
           0.000000     for  tau_reio
           0.000000     for  Omega0_scf
           0.000000     for  YHe
           0.000000     for  A_s
           0.000000     for  sigma8
           0.000000     for  100*theta_s
--> Not computing covariance matrix

I tried many different values of 'scf_parameters__2' and 'scf_parameters__5' with different permutations and combination and still I didn't get any valid output.

The .log file after the run look like this:

2023-02-08_200000__1.txt      Number of steps:200004     Steps accepted:810  acc = 0.004    min(-loglike) = 5850.91
2023-02-08_200000__2.txt      Number of steps:200000     Steps accepted:1119     acc = 0.0056   min(-loglike) = 5676.15
2023-02-08_200000__3.txt      Number of steps:200000     Steps accepted:1407     acc = 0.007    min(-loglike) = 5880.92
2023-02-08_200000__4.txt      Number of steps:200000     Steps accepted:1547     acc = 0.0077   min(-loglike) = 5941.33
2023-02-08_200000__5.txt      Number of steps:200000     Steps accepted:1181     acc = 0.0059   min(-loglike) = 5822.86
2023-02-08_200000__6.txt      Number of steps:200000     Steps accepted:1029     acc = 0.0051   min(-loglike) = 5510.55
2023-02-08_200000__7.txt      Number of steps:200000     Steps accepted:1579     acc = 0.0079   min(-loglike) = 5756.84
2023-02-08_200000__8.txt      Number of steps:200000     Steps accepted:1127     acc = 0.0056   min(-loglike) = 5667.24
2023-02-08_200000__9.txt      Number of steps:199998     Steps accepted:907  acc = 0.0045   min(-loglike) = 5517.73
2023-02-08_200000__10.txt     Number of steps:200000     Steps accepted:1583     acc = 0.0079   min(-loglike) = 5541.03
2023-02-08_200000__11.txt     Number of steps:200000     Steps accepted:1256     acc = 0.0063   min(-loglike) = 5630.97
2023-02-08_200000__12.txt     Number of steps:200000     Steps accepted:1189     acc = 0.0059   min(-loglike) = 5724.26
2023-02-08_200000__13.txt     Number of steps:196812     Steps accepted:1380     acc = 0.007    min(-loglike) = 5899.69
2023-02-08_200000__14.txt     Number of steps:200000     Steps accepted:927  acc = 0.0046   min(-loglike) = 5585.83
2023-02-08_200000__15.txt     Number of steps:200000     Steps accepted:1241     acc = 0.0062   min(-loglike) = 6442.02
2023-02-08_200000__16.txt     Number of steps:200000     Steps accepted:1617     acc = 0.0081   min(-loglike) = 6283.99
--> Total    number    of    steps: 3196815
--> Total number of accepted steps: 19900
--> Minimum of -logLike           : 5510.55

I have no clue how to solve the problem. It will be great if anyone can help me out on this.

Best, Amlan