California-Planet-Search / radvel

General Toolkit for Modeling Radial Velocity Data
http://radvel.readthedocs.io
MIT License
57 stars 52 forks source link

K2-24 tutorial hitting an error #316

Closed tronsgaard closed 4 years ago

tronsgaard commented 4 years ago

When I run the K2-24_Fitting+MCMC.ipynb tutorial, I hit an error after adding eccentricity priors in the last part:

res  = optimize.minimize(
    post.neglogprob_array, 
    post.get_vary_params(), 
    method='Powell',)

results in

---------------------------------------------------------------------------
AssertionError                            Traceback (most recent call last)
<ipython-input-15-e8c8bba411a1> in <module>
      2     post.neglogprob_array,
      3     post.get_vary_params(),
----> 4     method='Powell',)
      5 
      6 plot_results(like)

/usr/local/lib/python3.6/dist-packages/scipy/optimize/_minimize.py in minimize(fun, x0, args, method, jac, hess, hessp, bounds, constraints, tol, callback, options)
    598         return _minimize_neldermead(fun, x0, args, callback, **options)
    599     elif meth == 'powell':
--> 600         return _minimize_powell(fun, x0, args, callback, **options)
    601     elif meth == 'cg':
    602         return _minimize_cg(fun, x0, args, jac, callback, **options)

/usr/local/lib/python3.6/dist-packages/scipy/optimize/optimize.py in _minimize_powell(func, x0, args, callback, xtol, ftol, maxiter, maxfev, disp, direc, return_all, **unknown_options)
   2629         direc = asarray(direc, dtype=float)
   2630 
-> 2631     fval = squeeze(func(x))
   2632     x1 = x.copy()
   2633     iter = 0

/usr/local/lib/python3.6/dist-packages/scipy/optimize/optimize.py in function_wrapper(*wrapper_args)
    325     def function_wrapper(*wrapper_args):
    326         ncalls[0] += 1
--> 327         return function(*(wrapper_args + args))
    328 
    329     return ncalls, function_wrapper

/usr/local/lib/python3.6/dist-packages/radvel/likelihood.py in neglogprob_array(self, params_array)
    186 
    187     def neglogprob_array(self, params_array):
--> 188         return -self.logprob_array(params_array)
    189 
    190     def logprob_array(self, params_array):

/usr/local/lib/python3.6/dist-packages/radvel/posterior.py in logprob_array(self, param_values_array)
     61             float: log probability of the likelihood + priors
     62         """
---> 63         self.likelihood.set_vary_params(param_values_array)
     64         _logprob = self.logprob()
     65         # if not np.isfinite(_logprob):

/usr/local/lib/python3.6/dist-packages/radvel/likelihood.py in set_vary_params(self, param_values_array)
    148                 i += 1
    149             assert i == len(param_values_array), \
--> 150                 "Length of array must match number of varied parameters"
    151         except AttributeError:
    152             self.list_vary_params()

AssertionError: Length of array must match number of varied parameters

I'm using radvel version 1.4.0, Python 3.6.9

spencerhurt commented 4 years ago

I just went through the tutorial with version 1.4.0 and the optimization worked smoothly for me. It would probably help to see how you set up your posterior. If you want, you could download and post a pdf of your notebook and we can figure out what's going on.

tronsgaard commented 4 years ago

Thanks for looking into this. How strange. The only thing I modified was the path to the data. Happy to share the full notebook - please see http://tronsgaard.dk/files/K2-24_Fitting+MCMC.html

spencerhurt commented 4 years ago

I see what's happening! You're right, that's a bug on our end. When you're changing which parameters vary, the likelihood isn't recognizing the increase in varied params. In the meantime, it should work if you define a completely new likelihood object.

bjfultn commented 4 years ago

Thanks for looking into this @spencer3616!

bjfultn commented 4 years ago

This fix will be released in version 1.4.1

hposborn commented 8 months ago

Just an FYI for anyone who is still getting this error (NB, here's mine, python 3.9. Fresh conda environment. M2 Mac. radvel version 1.4.11)

IndexError                                Traceback (most recent call last)
Cell In[2], line 129
    126         print(simplepost)
    127 #         print(complxpost)
--> 129         simplemcmcout[star['Paper name']] = radvel.mcmc(simplepost, nwalkers=50, nrun=500)

File ~/miniconda3/envs/radvel/lib/python3.9/site-packages/radvel/mcmc.py:301, in mcmc(post, nwalkers, nrun, ensembles, checkinterval, minAfactor, maxArchange, burnAfactor, burnGR, maxGR, minTz, minsteps, minpercent, thin, serial, save, savename, proceed, proceedname, headless)
    298 # set up perturbation size
    300 pscales = []
--> 301 names = post.name_vary_params()
    302 for i,par in enumerate(post.vary_params):
    303     val = post.vector.vector[par][0]

File ~/miniconda3/envs/radvel/lib/python3.9/site-packages/radvel/likelihood.py:170, in Likelihood.name_vary_params(self)
    168 try:
    169     for i in self.vary_params:
--> 170         list.append(self.vector.names[i])
    171     return list
    172 except AttributeError:

IndexError: list index out of range

These are my parameter names using simplepost.list_vary_params:

                     value      vary
per1                        25.9049       True
tc1                         3003.44       True
secosw1                        0.01      False
sesinw1                        0.01      False
logk1                       1.55624       True
dvdt                     0.00273697       True
gamma_0                     18288.9       True
jit_0                       4.10976       True
tp1                          3000.2           
e1                           0.0002           
w1                         0.785398           
k1                          4.74095  

And therefore the indexes expected for varying parametrs is [0,1,4,5,6,7]. However, using simplepost.vary_params, I always get: array([0, 1, 4, 5, 7, 8]). So there is a "ghost variable" created which is not varied between parameters 5 and 6, and this shifts the indexes, resulting in the "list index out of range" and the code breaking when running the mcmc.

After much testing I realised this ghost parameter is 'curv' - i.e. the RV curvature. It is somehow included every time dvdt is initialised, though it does not appear in any of the parameter lists. In order to model a linearly-varying RV model and avoid this warning, it must be specified using params['curv'] = radvel.Parameter(value=0,vary=False). Maybe you can add an automatic initialisation at some point? Or remove this requirement...