search

kevin218 / Eureka

Eureka! is a data reduction and analysis pipeline intended for time-series observations with JWST.
https://eurekadocs.readthedocs.io/
MIT License
56 stars 44 forks source link

GP Model - ValueError: could not broadcast input array into shape #608

Closed sarthakpatel7 closed 5 months ago

sarthakpatel7 commented 5 months ago

Instrument

Light curve fitting (Stages 4-6)

What happened?

I am getting the following error when trying to use GP models with NIRSpec G395H data. GP model parameters are set with their default values and settings.

Error traceback output

Traceback (most recent call last):
  File "/mnt/data/spatel/DataAnalysis/g395h/input/nrs2/run_eureka.py", line 32, in <module>
    s5_meta = s5.fitlc(eventlabel, ecf_path=ecf_path) #, s4_meta=s4_meta)
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/eureka/S5_lightcurve_fitting/s5_fit.py", line 463, in fitlc
    meta, params = fit_channel(meta, time, flux, channel,
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/eureka/S5_lightcurve_fitting/s5_fit.py", line 859, in fit_channel
    lc_model.fit(model, meta, log, fitter='dynesty')
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/eureka/S5_lightcurve_fitting/lightcurve.py", line 179, in fit
    fit_model = self.fitter_func(self, model, meta, log, **kwargs)
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/eureka/S5_lightcurve_fitting/fitters.py", line 787, in dynestyfitter
    plots.plot_fit(lc, model, meta,
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/eureka/S5_lightcurve_fitting/plots_s5.py", line 54, in plot_fit
    model_gp = model.GPeval(model_noGP)
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/eureka/S5_lightcurve_fitting/models/Model.py", line 422, in GPeval
    flux = component.eval(fit, channel=channel,
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/eureka/S5_lightcurve_fitting/models/GPModel.py", line 186, in eval
    lcfinal = np.append(lcfinal, mu)
  File "<__array_function__ internals>", line 180, in append
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/numpy/lib/function_base.py", line 5390, in append
    values = ravel(values)
  File "<__array_function__ internals>", line 180, in ravel
  File "/home/spatel/anaconda3/envs/eureka/lib/python3.9/site-packages/numpy/core/fromnumeric.py", line 1859, in ravel
    return asanyarray(a).ravel(order=order)
ValueError: could not broadcast input array from shape (465,465) into shape (465,)

What operating system are you using?

Linux

What version of Python are you running?

Python 3.9

What Python packages do you have installed?

Name Version Build Channel

_libgcc_mutex 0.1 main
_openmp_mutex 5.1 1_gnu
alabaster 0.7.13 pypi_0 pypi asdf 2.15.0 pypi_0 pypi asdf-astropy 0.4.0 pypi_0 pypi asdf-coordinates-schemas 0.2.0 pypi_0 pypi asdf-standard 1.0.3 pypi_0 pypi asdf-transform-schemas 0.3.0 pypi_0 pypi asdf-unit-schemas 0.1.0 pypi_0 pypi asdf-wcs-schemas 0.1.1 pypi_0 pypi asteval 0.9.29 pypi_0 pypi astraeus 0.3 pypi_0 pypi astropy 5.2.2 pypi_0 pypi astropy-healpix 0.7 pypi_0 pypi astroquery 0.4.6 pypi_0 pypi astroscrappy 1.1.0 pypi_0 pypi asttokens 2.2.1 pypi_0 pypi attrs 23.1.0 pypi_0 pypi babel 2.12.1 pypi_0 pypi backcall 0.2.0 pypi_0 pypi batman-package 2.4.9 pypi_0 pypi bayesicfitting 3.1.1 pypi_0 pypi beautifulsoup4 4.12.2 pypi_0 pypi bokeh 2.4.3 pypi_0 pypi bottleneck 1.3.7 pypi_0 pypi ca-certificates 2023.01.10 h06a4308_0
ccdproc 2.4.0 pypi_0 pypi celerite 0.4.2 pypi_0 pypi certifi 2022.5.18.1 pypi_0 pypi cffi 1.15.1 pypi_0 pypi cftime 1.6.2 pypi_0 pypi charset-normalizer 3.1.0 pypi_0 pypi cloudpickle 2.2.1 pypi_0 pypi contourpy 1.0.7 pypi_0 pypi corner 2.2.2 pypi_0 pypi crds 11.17.0 pypi_0 pypi cryptography 40.0.2 pypi_0 pypi cycler 0.11.0 pypi_0 pypi cython 0.29.34 pypi_0 pypi dask 2022.6.0 pypi_0 pypi decorator 5.1.1 pypi_0 pypi docutils 0.20 pypi_0 pypi drizzle 1.13.7 pypi_0 pypi dynesty 2.1.1 pypi_0 pypi emcee 3.1.4 pypi_0 pypi eureka 0.11.dev35+gfb357182 pypi_0 pypi exceptiongroup 1.1.1 pypi_0 pypi executing 1.2.0 pypi_0 pypi exotic-ld 3.0.0 pypi_0 pypi filelock 3.12.0 pypi_0 pypi fonttools 4.39.4 pypi_0 pypi fsspec 2023.5.0 pypi_0 pypi future 0.18.3 pypi_0 pypi george 0.4.0 pypi_0 pypi gwcs 0.18.3 pypi_0 pypi h5netcdf 1.1.0 pypi_0 pypi h5py 3.1.0 pypi_0 pypi html5lib 1.1 pypi_0 pypi idna 3.4 pypi_0 pypi imageio 2.28.1 pypi_0 pypi imagesize 1.4.1 pypi_0 pypi importlib-metadata 6.6.0 pypi_0 pypi importlib-resources 5.12.0 pypi_0 pypi iniconfig 2.0.0 pypi_0 pypi ipython 8.13.2 pypi_0 pypi jaraco-classes 3.2.3 pypi_0 pypi jedi 0.18.2 pypi_0 pypi jeepney 0.8.0 pypi_0 pypi jinja2 3.1.2 pypi_0 pypi jmespath 1.0.1 pypi_0 pypi jsonschema 4.9.1 pypi_0 pypi jwst 1.11.4 pypi_0 pypi keyring 23.13.1 pypi_0 pypi kiwisolver 1.4.4 pypi_0 pypi lazy-loader 0.2 pypi_0 pypi ld_impl_linux-64 2.38 h1181459_1
libffi 3.3 he6710b0_2
libgcc-ng 11.2.0 h1234567_1
libgomp 11.2.0 h1234567_1
libstdcxx-ng 11.2.0 h1234567_1
lmfit 1.2.1 pypi_0 pypi locket 1.0.0 pypi_0 pypi markupsafe 2.1.2 pypi_0 pypi matplotlib 3.7.1 pypi_0 pypi matplotlib-inline 0.1.6 pypi_0 pypi more-itertools 9.1.0 pypi_0 pypi ncurses 6.4 h6a678d5_0
netcdf4 1.6.3 pypi_0 pypi networkx 3.1 pypi_0 pypi numpy 1.22.0 pypi_0 pypi numpydoc 1.5.0 pypi_0 pypi opencv-python-headless 4.8.1.78 pypi_0 pypi openssl 1.1.1t h7f8727e_0
packaging 23.1 pypi_0 pypi pandas 2.0.1 pypi_0 pypi parsley 1.3 pypi_0 pypi parso 0.8.3 pypi_0 pypi partd 1.4.0 pypi_0 pypi pexpect 4.8.0 pypi_0 pypi photutils 1.7.0 pypi_0 pypi pickleshare 0.7.5 pypi_0 pypi pillow 9.5.0 pypi_0 pypi pip 23.0.1 py39h06a4308_0
pluggy 1.0.0 pypi_0 pypi poppy 1.1.0 pypi_0 pypi prompt-toolkit 3.0.38 pypi_0 pypi psutil 5.9.5 pypi_0 pypi ptyprocess 0.7.0 pypi_0 pypi pure-eval 0.2.2 pypi_0 pypi pycparser 2.21 pypi_0 pypi pyerfa 2.0.0.3 pypi_0 pypi pygments 2.15.1 pypi_0 pypi pyparsing 3.0.9 pypi_0 pypi pyrsistent 0.19.3 pypi_0 pypi pysynphot 2.0.0 pypi_0 pypi pytest 7.3.1 pypi_0 pypi python 3.9.7 h12debd9_1
python-dateutil 2.8.2 pypi_0 pypi pytz 2023.3 pypi_0 pypi pyvo 1.4.1 pypi_0 pypi pywavelets 1.4.1 pypi_0 pypi pyyaml 6.0 pypi_0 pypi readline 8.2 h5eee18b_0
reproject 0.10.0 pypi_0 pypi requests 2.30.0 pypi_0 pypi scikit-image 0.20.0 pypi_0 pypi scipy 1.9.1 pypi_0 pypi secretstorage 3.3.3 pypi_0 pypi semantic-version 2.10.0 pypi_0 pypi setuptools 66.0.0 py39h06a4308_0
setuptools-scm 7.1.0 pypi_0 pypi six 1.16.0 pypi_0 pypi snowballstemmer 2.2.0 pypi_0 pypi soupsieve 2.4.1 pypi_0 pypi spherical-geometry 1.2.23 pypi_0 pypi sphinx 7.0.1 pypi_0 pypi sphinxcontrib-applehelp 1.0.4 pypi_0 pypi sphinxcontrib-devhelp 1.0.2 pypi_0 pypi sphinxcontrib-htmlhelp 2.0.1 pypi_0 pypi sphinxcontrib-jsmath 1.0.1 pypi_0 pypi sphinxcontrib-qthelp 1.0.3 pypi_0 pypi sphinxcontrib-serializinghtml 1.1.5 pypi_0 pypi sqlite 3.41.2 h5eee18b_0
stack-data 0.6.2 pypi_0 pypi stcal 1.4.4 pypi_0 pypi stdatamodels 1.7.2 pypi_0 pypi stpipe 0.5.0 pypi_0 pypi stsci-image 2.3.5 pypi_0 pypi stsci-imagestats 1.6.3 pypi_0 pypi stsci-stimage 0.2.6 pypi_0 pypi svo-filters 0.4.4 pypi_0 pypi tifffile 2023.4.12 pypi_0 pypi tk 8.6.12 h1ccaba5_0
tomli 2.0.1 pypi_0 pypi toolz 0.12.0 pypi_0 pypi tornado 6.3.2 pypi_0 pypi tqdm 4.65.0 pypi_0 pypi traitlets 5.9.0 pypi_0 pypi tweakwcs 0.8.2 pypi_0 pypi typing-extensions 4.5.0 pypi_0 pypi tzdata 2023.3 pypi_0 pypi uncertainties 3.1.7 pypi_0 pypi urllib3 2.0.2 pypi_0 pypi wcwidth 0.2.6 pypi_0 pypi webencodings 0.5.1 pypi_0 pypi wheel 0.38.4 py39h06a4308_0
wiimatch 0.3.1 pypi_0 pypi xarray 2023.4.2 pypi_0 pypi xz 5.4.2 h5eee18b_0
zipp 3.15.0 pypi_0 pypi zlib 1.2.13 h5eee18b_0

Code of Conduct

taylorbell57 commented 5 months ago

Strange, I use the GP model all the time. Would you mind copy-pasting the contents of your Stage 5 ECF and EPF so that I have more information with which to troubleshoot?

sarthakpatel7 commented 5 months ago 
# Eureka! Control File for Stage 5: Lightcurve Fitting

# Stage 5 Documentation: https://eurekadocs.readthedocs.io/en/latest/ecf.html#stage-5

ncpu            50 # The number of CPU threads to use when running emcee or dynesty in parallel

allapers        False                   # Run S5 on all of the apertures considered in S4? Otherwise will use newest output in the inputdir
rescale_err     False                   # Rescale uncertainties to have reduced chi-squared of unity
fit_par         ./S5_fit_par_wasp39b.epf   # What fitting epf do you want to use?
verbose         True                    # If True, more details will be printed about steps
fit_method      [dynesty]               #options are: lsq, emcee, dynesty (can list multiple types separated by commas)
run_myfuncs     [batman_tr,polynomial,step]  # options are: batman_tr, batman_ecl, sinusoid_pc, expramp, polynomial, step, xpos, ypos, xwidth, ywidth, and GP (can list multiple types separated by commas)

# Manual clipping in time
manual_clip     None    # A list of lists specifying the start and end integration numbers for manual removal.

# Limb darkening controls
# IMPORTANT: limb-darkening coefficients are not automatically fixed then, change to 'fixed' in .epf file whether they should be fixed or fitted!
use_generate_ld  exotic-ld  # use the generated limb-darkening coefficients from Stage 4? Options: exotic-ld, None. For exotic-ld, the limb-darkening laws available are linear, quadratic, 3-parameter and 4-parameter non-linear.
ld_file          None  # Fully qualified path to the location of a limb darkening file that you want to use
ld_file_white    None  # Fully qualified path to the location of a limb darkening file that you want to use for the white-light light curve (required if ld_file is not None and any EPF parameters are set to white_free or white_fixed).

# General fitter
old_fitparams   None # filename relative to topdir that points to a fitparams csv to resume where you left off (set to None to start from scratch)

#lsq
lsq_method      'Powell' # The scipy.optimize.minimize optimization method to use
lsq_tol         1e-6    # The tolerance for the scipy.optimize.minimize optimization method
lsq_maxiter     None    # Maximum number of iterations to perform. Depending on the method each iteration may use several function evaluations. Set to None to use the default value

#mcmc
old_chain       None    # Output folder relative to topdir that contains an old emcee chain to resume where you left off (set to None to start from scratch)
lsq_first       True    # Initialize with an initial lsq call (can help shorten burn-in, but turn off if lsq fails). Only used if old_chain is None
run_nsteps      1000
run_nwalkers    200
run_nburn       500     # How many of run_nsteps should be discarded as burn-in steps

#dynesty
run_nlive       512    # Must be > ndim * (ndim + 1) // 2
run_bound       'multi'
run_sample      'auto'
run_tol         0.1

#GP inputs
kernel_inputs   ['time'] #options: time
kernel_class    ['Matern32'] #options: ExpSquared, Matern32, Exp, RationalQuadratic for george, Matern32 for celerite (sums of kernels possible for george separated by commas)
GP_package      'celerite' #options: george, celerite
useHODLR        False

# Plotting controls
interp          False   # Should astrophysical model be interpolated (useful for uneven sampling like that from HST)

# Diagnostics
isplots_S5      5       # Generate few (1), some (3), or many (5) figures (Options: 1 - 5)
nbin_plot       100     # The number of bins that should be used for figures 5104 and 5304. Defaults to 100.
testing_S5      False   # Boolean, set True to only use the first spectral channel
testing_model   False   # Boolean, set True to only inject a model source of systematics
hide_plots      True   # If True, plots will automatically be closed rather than popping up

# Project directory
topdir          /mnt/data/spatel/DataAnalysis/g395h/output/

# Directories relative to topdir
inputdir        Stage4/nrs2/S4_2024-01-25_wasp39b_run2   # The folder containing the outputs from Eureka!'s S4 pipeline (will be overwritten if calling S4 and S5 sequentially)
outputdir       Stage5/nrs2
sarthakpatel7 commented 5 months ago 
# Stage 5 Fit Parameters Documentation: https://eurekadocs.readthedocs.io/en/latest/ecf.html#stage-5-fit-parameters

# Name       Value         Free?          PriorPar1    PriorPar2    PriorType
# "Free?" can be free, fixed, white_free, white_fixed, shared, or independent
# PriorType can be U (Uniform), LU (Log Uniform), or N (Normal).
# If U/LU, PriorPar1 and PriorPar2 represent upper and lower limits of the parameter/log(the parameter).
# If N, PriorPar1 is the mean and PriorPar2 is the standard deviation of a Gaussian prior.
#-------------------------------------------------------------------------------------------------------
#
# ------------------
# ** Transit/eclipse parameters **
# ------------------
rp           0.148          'free'         0.12         0.16          U
#fp           0.008          'free'         0           0.5          U
# ----------------------
# ** Phase curve parameters **
# ----------------------
#AmpCos1      0.4           'free'         0            1            U
#AmpSin1      0.01          'free'         -1           1            U
#AmpCos2      0.01          'free'         -1           1            U
#AmpSin2      0.01          'free'         -1           1            U
# ------------------
# ** Orbital parameters **
# ------------------
per          4.055259       'fixed'         4     4.1         U
t0           57422.840      'free'         57000    61000         U
time_offset  0              'independent'
inc          87.78          'free'         83        90         U
a            11.485           'free'         11.485         0.015            N  # An incorrect a/R* seems to have been used in the simulations
ecc          0.0            'fixed'        0            1            U
w            90.            'fixed'        0            180          U
# -------------------------
# ** Limb darkening parameters **
# Choose limb_dark from ['uniform', 'linear', 'quadratic', 'kipping2013', 'squareroot', 'logarithmic', 'exponential','3-parameter', '4-parameter']
# When using generated limb-darkening coefficients from exotic-ld choose from ['linear', 'quadratic', '3-parameter', '4-parameter']
# -------------------------
limb_dark    'quadratic'  'independent'
u1           0.59           'free'         0            1            U
u2           0.53           'fixed'         0            1            U
# --------------------
# ** Systematic variables **
# Polynomial model variables (c0--c9 for 0th--3rd order polynomials in time); Fitting at least c0 is very strongly recommended!
# Exponential ramp model variables (r0--r2 for one exponential ramp, r3--r5 for a second exponential ramp)
# GP model parameters (A, m for the first kernel; A1, m1 for the second kernel; etc.) in log scale
A            -7             'free'         -15          0            U
m            -4             'free'         -6           0            U
# Step-function model variables (step# and steptime# for step-function model #; e.g. step0 and steptime0)
# Drift model variables (xpos, ypos, xwidth, ywidth)
# --------------------
c0           1.015          'free'         0.8            1.2         U
c1           0              'free'         0            0.01         N
# -----------
# ** White noise **
# Use scatter_mult to fit a multiplier to the expected noise level from Stage 3 (recommended)
# Use scatter_ppm to fit the noise level in ppm
# -----------
scatter_mult 1.1            'free'         1.1            10          N
kevin218 commented 5 months ago

Under run_myfuncs, it looks like you're using 'step' and not 'GP'. I don't know how you're even getting this error without 'GP'.

Unrelated, you may want to fix both LD parameters when using exotic-ld.

On Tue, Jan 30, 2024, 5:27 AM sarthakpatel7 @.***> wrote:

Stage 5 Fit Parameters Documentation: https://eurekadocs.readthedocs.io/en/latest/ecf.html#stage-5-fit-parameters Name Value Free? PriorPar1 PriorPar2 PriorType "Free?" can be free, fixed, white_free, white_fixed, shared, or independent PriorType can be U (Uniform), LU (Log Uniform), or N (Normal). If U/LU, PriorPar1 and PriorPar2 represent upper and lower limits of the parameter/log(the parameter). If N, PriorPar1 is the mean and PriorPar2 is the standard deviation of a Gaussian prior.

-------------------------------------------------------------------------------------------------------

------------------ Transit/eclipse parameters ------------------

rp 0.148 'free' 0.12 0.16 U

fp 0.008 'free' 0 0.5 U

---------------------- Phase curve parameters ----------------------

AmpCos1 0.4 'free' 0 1 U

AmpSin1 0.01 'free' -1 1 U

AmpCos2 0.01 'free' -1 1 U

AmpSin2 0.01 'free' -1 1 U

------------------ Orbital parameters ------------------

per 4.055259 'fixed' 4 4.1 U t0 57422.840 'free' 57000 61000 U time_offset 0 'independent' inc 87.78 'free' 83 90 U a 11.485 'free' 11.485 0.015 N # An incorrect a/R* seems to have been used in the simulations ecc 0.0 'fixed' 0 1 U w 90. 'fixed' 0 180 U ------------------------- Limb darkening parameters Choose limb_dark from ['uniform', 'linear', 'quadratic', 'kipping2013', 'squareroot', 'logarithmic', 'exponential','3-parameter', '4-parameter'] When using generated limb-darkening coefficients from exotic-ld choose from ['linear', 'quadratic', '3-parameter', '4-parameter']

limb_dark 'quadratic' 'independent' u1 0.59 'free' 0 1 U u2 0.53 'fixed' 0 1 U -------------------- Systematic variables Polynomial model variables (c0--c9 for 0th--3rd order polynomials in time); Fitting at least c0 is very strongly recommended! Exponential ramp model variables (r0--r2 for one exponential ramp, r3--r5 for a second exponential ramp) GP model parameters (A, m for the first kernel; A1, m1 for the second kernel; etc.) in log scale

A -7 'free' -15 0 U m -4 'free' -6 0 U Step-function model variables (step# and steptime# for step-function model

; e.g. step0 and steptime0) Drift model variables (xpos, ypos, xwidth,

ywidth) --------------------

c0 1.015 'free' 0.8 1.2 U c1 0 'free' 0 0.01 N ----------- White noise Use scatter_mult to fit a multiplier to the expected noise level from Stage 3 (recommended) Use scatter_ppm to fit the noise level in ppm -----------

scatter_mult 1.1 'free' 1.1 10 N

— Reply to this email directly, view it on GitHub https://github.com/kevin218/Eureka/issues/608#issuecomment-1916535769, or unsubscribe https://github.com/notifications/unsubscribe-auth/AFC2C7OEUYCFP3SWWDSTMK3YRDDJRAVCNFSM6AAAAABCOQMRB2VHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMYTSMJWGUZTKNZWHE . You are receiving this because you are subscribed to this thread.Message ID: @.***>

sarthakpatel7 commented 5 months ago

It seems I copied the wrong ecf file. It was 'GP' instead of 'step' in under run_myfuncs. Apologies for the mistake.

taylorbell57 commented 5 months ago

Can you confirm that you only get this bug when you're fitting for the GP model and not otherwise?

taylorbell57 commented 5 months ago

Aha, I was just able to reproduce this error myself - I should have a patch ready in the next hour or so