MNGuenther / allesfitter

allesfitter is a convenient wrapper around the packages ellc (light curve and RV models), dynesty (static and dynamic nested sampling) emcee (Markov Chain Monte Carlo sampling) and celerite (Gaussian Process models).
MIT License
60 stars 36 forks source link

Weird depth and rr fit values vs params_star #22

Closed martindevora closed 3 years ago

martindevora commented 3 years ago

Hi. First of all I don't want to state this issue as a bug before gathering info from you. Let me expose what I'm trying: I'm running allesfitter with: allesfitter.ns_fit(allesfit_dir) allesfitter.ns_output(allesfit_dir)

I'm providing settings.csv, params.csv with this content:

#companion b astrophysical params,,,,,
SOI_5_rr,0.026,1,trunc_normal 0 1 0.026 0.05,$R_b / R_\star$,
SOI_5_rsuma,0.11011430055938672,1,trunc_normal 0 1 0.11011430055938672 0.1,$(R_\star + R_b) / a_b$,
SOI_5_cosi,0.0,1,uniform 0.0 0.06,$\cos{i_b}$,
SOI_5_epoch,1362.8847530987089,1,trunc_normal -1000000000000.0 1000000000000.0 1362.8847530987089 0.05,$T_{0;b}$,$\mathrm{BJD}$
SOI_5_period,13.189658223908912,1,trunc_normal -1000000000000.0 1000000000000.0 13.189658223908912 0.05,$P_b$,$\mathrm{d}$
SOI_5_f_c,0.0,0,trunc_normal -1 1 0.0 0.0,$\sqrt{e_b} \cos{\omega_b}$,
SOI_5_f_s,0.0,0,trunc_normal -1 1 0.0 0.0,$\sqrt{e_b} \sin{\omega_b}$,

#dilution per instrument,,,,,
dil_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$D_\mathrm{0; lc}$,
#limb darkening coefficients per instrument,,,,,
host_ldc_q1_lc,0.1998,0,uniform 0.0 1.0,$q_{1; \mathrm{lc}}$,
host_ldc_q2_lc,0.5387,0,uniform 0.0 1.0,$q_{2; \mathrm{lc}}$,
#surface brightness per instrument and companion,,,,,
SOI_5_sbratio_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$J_{b; \mathrm{lc}}$,
#albedo per instrument and companion,,,,,
host_geom_albedo_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$A_{\mathrm{geom}; host; \mathrm{lc}}$,
SOI_5_geom_albedo_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$A_{\mathrm{geom}; b; \mathrm{lc}}$,
#gravity darkening per instrument and companion,,,,,
host_gdc_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$Grav. dark._{b; \mathrm{lc}}$,
#spots per instrument and companion,,,,,
#errors per instrument,
ln_err_flux_lc,-7.0,1,uniform -15.0 0.0,$\log{\sigma_\mathrm{lc}}$,$\log{ \mathrm{rel. flux.} }$
#baseline per instrument,
baseline_gp_matern32_lnsigma_flux_lc,0.0,1,uniform -15.0 15.0,$\mathrm{gp: \ln{\sigma} (lc)}$,
baseline_gp_matern32_lnrho_flux_lc,0.0,1,uniform -15.0 15.0,$\mathrm{gp: \ln{\rho} (lc)}$,

As you can see, I'm providing a period of 13.19d and a rr of 0.026.

and params_star.csv with content:

165.691476612015,-16.4061904736575,1.133,0.0645,0.0645,0.992,0.082,0.082,5540.0,220,220,0.1998,0.5387

As you can see, I'm providing a R_star of 1.133 with little uncertainty.

Given these priors, I obtain the next posteriors:

Companion radius SOI_5; $R_\mathrm{SOI_5}$
($\mathrm{R_{\oplus}}$),12.997667627520098,3.6315454554099897,4.360371356068896,derived
SOI_5_rr,0.10548749875981736,0.02878261480797195,0.033790848181602606,$R_b / R_\star$,
Transit depth (dil.) SOI_5; $\delta_\mathrm{tr; dil; SOI_5; lc}$ (ppt),0.0012701552787879322,0.0001110527739457412,0.00011386618070774457,derived

Which seems to be unlikely with the star params I provided because: With the obtained depth of 1 ppt (0.00127) the R_b/R_star value would be aproximately 0.04 (given sqrt(depth) = R_b/R_star). Am I missing anything? I know I'm not restricting the uncertainties from params.csv but I'd expect that the star parameters that I'm providing as input would act as a hard limit.

Kind regards and again, thank you in advance for your attention.

MNGuenther commented 3 years ago

Hi @martindevora, thanks for reaching out!

(1) I'm a bit confused about your priors here (also see attached figure): SOI_5_rr,0.026,1,trunc_normal 0 1 0.07 0.05,$Rb / R\star$, Yes, you give a starting guess of 0.026 for the fit, but your prior enforces a mean of 0.07 with std of 0.05. How come?

(2) Generally, one cannot use normal priors unless that info comes from somewhere. For example, if you have fitted TESS data before, you may then use the TESS posteriors as priors in your fit of SOI data. But if you have not done that before, the SOI priors should all be uniform.

(3) The transit depth provided by allesfitter is not Rb^2/Rstar^2. Instead, it is the actually measured transit depth, including all effects of grazing-transits (cos i), limb darkening and dilution.

Without having seen the full posteriors or plots, my first guess is that your prior on rr enforced a solution where rr is large. This gets compensated by a large cos i (grazing transit) in the fit. If I am right, you should see a large value for cos i and a strong correlation between rr and cos i in your corner plot.

Feel free to email me the full folder (or at least all plots and tables) and I can have a closer look.

Cheers Max image

MNGuenther commented 3 years ago

Btw (and this is my fault for not polishing the GUI more): you can delete the "useless" (not-used) rows from your params.csv file, which will make it a bit cleaner. Such as these:

#dilution per instrument,,,,,
dil_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$D_\mathrm{0; lc}$,
[...]
#surface brightness per instrument and companion,,,,,
SOI_5_sbratio_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$J_{b; \mathrm{lc}}$,
#albedo per instrument and companion,,,,,
host_geom_albedo_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$A_{\mathrm{geom}; host; \mathrm{lc}}$,
SOI_5_geom_albedo_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$A_{\mathrm{geom}; b; \mathrm{lc}}$,
#gravity darkening per instrument and companion,,,,,
host_gdc_lc,0.0,0,trunc_normal 0 1 0.0 0.0,$Grav. dark._{b; \mathrm{lc}}$,
#spots per instrument and companion,,,,
martindevora commented 3 years ago

Sorry. I copied one file content before the final one and something got mixed in between, but I edited the post 5 minutes after creating it and I see the mean has the same value. I don't know why you are still seeing the original post content, sorry for that.

SOI_5_rr,0.026,1,trunc_normal 0 1 0.026 0.05,$R_b / R_\star$

I'm already using the priors obtained from a TLS execution. However, as the depth comes without an error, I'm giving the R_planet/R_star ratio a large uncertainty.

I don't know what I did with the corner plot for these parameters, but I can share another one, which is very similar, where I varied a little bit the priors and the result is pretty much the same: ns_corner.pdf

I guess that the correlation can be appreciated between the cos(i) and (R*+Rp / a)?

As far as I understand (sorry I'm just a beginner with transit modelling) the model is guessing a much larger radius because it is interpreting the transit as a grazing one, which would have made the depth to be shallower. That'd make sense even though it is not the result that I'm expecting from the visual inspection.

Thank you for pointing to that.

MNGuenther commented 3 years ago

I'm already using the priors obtained from a TLS execution.

That raises alarm bells! One cannot re-use the same data to get priors and then do a fit! In Bayesian statistics, you can only use the data once. For example, you can use TESS data to get a prior, and then use that prior to fit SOI data (I don't know the abbreviation SOI, I guess it's a telescope?). But one can never-ever run BLS/TLS on SOI data and then use the outcome as a prior to fit the same SOI data!

You can (and should) use the TLS outcome for your initial guess, but you must leave the priors uniform everywhere.

Could you re-run the fit after those fixes and send me the full directory and results?

martindevora commented 3 years ago

Sorry. SOI is a label that we gave to some TESS data.

We are running TLS on TESS data to get an approximation for period, epoch and depth (rp/rstar). The bayesian fit is then executed against the same TESS data given those restricted priors, so we can get refined values of them and of course, the posteriors. The final goal is to be able to plan observations based on the ephemeris obtained from the final priors.

So you are telling me to always use "uniform" for all the priors? If so, I'll gladly re-run them and share the results to you.

Thanks and regards.

MNGuenther commented 3 years ago

We are running TLS on TESS data to get an approximation for period, epoch and depth (rp/rstar). The bayesian fit is then executed against the same TESS data given those restricted priors, so we can get refined values of them and of course, the posteriors.

Exactly, that is something one can not do with Bayesian statistics. Since you use the same data set (TESS), the correct way to do it is: (1) run TLS to find the signal and get an initial idea of period, epoch and depth; (2) use that as initial guess in your MCMC fit; (3) use uniform priors.

You can only use normal priors if you have a reason for it that comes from some other data. For example, if you had Kepler results for that same TESS target but you don't want to re-analyzie Kepler, you could use the Kepler posteriors as priors for your TESS analysis.

martindevora commented 3 years ago

Now I'm understanding you better. Then, I'll change my priors distributions and keep using the initial guess.

Thanks for pointing at the root of the issue and being patient with your answers.

Regards.

MNGuenther commented 3 years ago

My pleasure and glad I could help :) I don't think this will give you a very different radius, so if you are still puzzled by your new results (after re-running), feel free to email me the folder and I'll have a look.