timothydmorton / isochrones

Pythonic stellar model grid access; easy MCMC fitting of stellar properties
http://isochrones.readthedocs.org
MIT License
123 stars 63 forks source link

Age prior not being respected #174

Closed warrickball closed 2 years ago

warrickball commented 2 years ago

I'm dealing with some low-mass stars and getting skewed results because of the default age prior that bounds log₁₀(age) between 5 and 10.15 (which is about the age of the Universe). The MIST isochrones go up to log₁₀(age)=10.3, which I've tested by generating the isochrones so I tried changing the age prior with

from isochrones import get_ichrone, SingleStarModel, priors

mist = get_ichrone('mist', bands=['G', 'BP', 'RP'])

params = {'Teff': (4790, 77), 'feh': (-0.3, 0.1), 'logg': (3.15, 0.06),
          'parallax': (10.1508, 0.5)}  # mas                                                                                                                                                                                                                  
mod = SingleStarModel(mist, **params)
mod._priors['AV'] = priors.PowerLawPrior(-3, bounds=(0.001, 1.0)) # like a δ-function at Av=0
mod._priors['mass'] = priors.FlatPrior((0.5, 3.5))
mod._priors['feh'] = priors.FlatPrior((-2.0, 1.0))
mod._priors['distance'] = priors.DistancePrior(max_distance=1000)
mod._priors['age'] = priors.FlatLogPrior((7,10.3))
mod.fit()

fig = mod.corner_physical()
fig.savefig('physical.png')
fig = mod.corner_observed()
fig.savefig('observed.png')
print(mod.derived_samples.to_records()['age'].max())

but my fits all still seem to be restricted to log₁₀(age)<10.15. E.g., this script's terminal output from the last line is always just under 10.15, rather than closer to 10.3. Corner plot of parameters. How should I be extending the age range? (I don't mind editing the package source itself.)

warrickball commented 2 years ago

I found the bounds hardcoded in mist/models.py:

class MISTModelGrid(StellarModelGrid):
    ...
    bounds = (("age", (5, 10.15)), ("feh", (-4, 0.5)), ("eep", (0, 1710)), ("mass", (0.1, 300)))
    ...

I've raised to the upper bound on age to 10.27697, which seems to have been the maximum isochrone I could generate with the mod.ic.isochrone function, and I now get older (and less massive) models.