From 7ce3cdb10fb3d16f4758b2430f1e496ee37fd442 Mon Sep 17 00:00:00 2001 From: Richard West Date: Mon, 4 May 2026 09:25:26 -0400 Subject: [PATCH 1/9] Update Arkane prettifier for Constant AST nodes Python 3.11 parses literal strings and numbers as ast.Constant nodes and warns when NodeVisitor falls back to visit_Str or visit_Num. Updating the prettifier to recognize Constant nodes removes those deprecation warnings while keeping the existing pretty-printing logic localized to Arkane output formatting. This is a good fit because it follows the modern AST API directly instead of filtering warnings or depending on deprecated compatibility dispatch. The change is intentionally narrow: it only affects literal classification and rendering inside prettify(). Co-authored-by: Codex with GPT-5.5 --- arkane/output.py | 44 ++++++++++++++++++++++++++------------------ 1 file changed, 26 insertions(+), 18 deletions(-) diff --git a/arkane/output.py b/arkane/output.py index 7657834c004..1836229c765 100644 --- a/arkane/output.py +++ b/arkane/output.py @@ -56,6 +56,18 @@ def __init__(self, level=0, indent=4): self.level = level self.indent = indent + @staticmethod + def _is_string(node): + return isinstance(node, ast.Constant) and isinstance(node.value, str) + + @staticmethod + def _is_number(node): + return isinstance(node, ast.Constant) and isinstance(node.value, (int, float, complex)) and not isinstance(node.value, bool) + + @classmethod + def _is_string_or_number(cls, node): + return cls._is_string(node) or cls._is_number(node) + def visit_Call(self, node): """ Return a pretty representation of the class or function call represented by `node`. @@ -84,7 +96,7 @@ def visit_List(self, node): """ Return a pretty representation of the list represented by `node`. """ - if any([not isinstance(e, (ast.Str, ast.Num, ast.UnaryOp)) for e in node.elts]): + if any([not (self._is_string_or_number(e) or isinstance(e, ast.UnaryOp)) for e in node.elts]): # Split elements onto multiple lines result = '[\n' self.level += 1 @@ -105,11 +117,11 @@ def visit_Tuple(self, node): """ # If the tuple represents a quantity, keep it on one line is_quantity = True - if len(node.elts) == 0 or not isinstance(node.elts[0], (ast.Num, ast.List)) or ( + if len(node.elts) == 0 or not (self._is_number(node.elts[0]) or isinstance(node.elts[0], ast.List)) or ( isinstance(node.elts[0], ast.List) and - any([not isinstance(e, (ast.Num, ast.UnaryOp)) for e in node.elts[0].elts])): + any([not (self._is_number(e) or isinstance(e, ast.UnaryOp)) for e in node.elts[0].elts])): is_quantity = False - elif len(node.elts) < 2 or not isinstance(node.elts[1], ast.Str): + elif len(node.elts) < 2 or not self._is_string(node.elts[1]): is_quantity = False if not is_quantity: @@ -131,8 +143,8 @@ def visit_Dict(self, node): """ Return a pretty representation of the dict represented by `node`. """ - if (any([not isinstance(e, (ast.Str, ast.Num)) for e in node.keys]) - or any([not isinstance(e, (ast.Str, ast.Num)) for e in node.values])): + if (any([not self._is_string_or_number(e) for e in node.keys]) + or any([not self._is_string_or_number(e) for e in node.values])): # Split elements onto multiple lines result = '{\n' self.level += 1 @@ -149,20 +161,16 @@ def visit_Dict(self, node): self.string = result return result - def visit_Str(self, node): - """ - Return a pretty representation of the string represented by `node`. - """ - result = repr(node.s) - self.string = result - return result - - def visit_Num(self, node): + def visit_Constant(self, node): """ - Return a pretty representation of the number represented by `node`. + Return a pretty representation of the constant represented by `node`. """ - result = '{0:g}'.format(node.n) - # result = repr(node.n) + if isinstance(node.value, str): + result = repr(node.value) + elif isinstance(node.value, (int, float, complex)) and not isinstance(node.value, bool): + result = '{0:g}'.format(node.value) + else: + return None self.string = result return result From 6a88539090bc2fe7579aefbded8e0dac125b8fa6 Mon Sep 17 00:00:00 2001 From: Richard West Date: Mon, 4 May 2026 09:25:36 -0400 Subject: [PATCH 2/9] Handle underdetermined BAC confidence intervals explicitly BAC cross-validation can fit folds with no residual degrees of freedom, so the variance estimate is mathematically undefined. Letting NumPy perform the division produced runtime warnings for an expected edge case and made warning-clean test runs noisy. The implementation makes the undefined case explicit by returning nan confidence interval values when n - p is not positive, and by masking invalid covariance-to-correlation divisions. This preserves the numerical meaning of the result while avoiding warning-driven control flow. Consideration: callers that inspect confidence intervals from underdetermined folds should continue to treat nan values as unavailable uncertainty estimates. Co-authored-by: Codex with GPT-5.5 Co-authored-by: Jackson Burns <33505528+JacksonBurns@users.noreply.github.com> --- arkane/encorr/bac.py | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/arkane/encorr/bac.py b/arkane/encorr/bac.py index c511a9f47b4..da1cc5d9075 100644 --- a/arkane/encorr/bac.py +++ b/arkane/encorr/bac.py @@ -1092,10 +1092,13 @@ def get_confidence_intervals(x: np.ndarray, weights = np.eye(n) e = y - ypred # Residuals - sigma2 = e.T @ weights @ e / (n - p) # MSE + dof = n - p + if dof <= 0: + return np.nan, np.nan + sigma2 = e.T @ weights @ e / dof # MSE cov = sigma2 * np.linalg.inv(x.T @ weights @ x) # covariance matrix se = np.sqrt(np.diag(cov)) # standard error - tdist = distributions.t.ppf(1 - alpha / 2, n - p) # student-t + tdist = distributions.t.ppf(1 - alpha / 2, dof) # student-t ci = tdist * se # confidence interval half-width return ci, cov @@ -1103,6 +1106,7 @@ def get_confidence_intervals(x: np.ndarray, def _covariance_to_correlation(cov: np.ndarray) -> np.ndarray: """Convert (unscaled) covariance matrix to correlation matrix""" v = np.sqrt(np.diag(cov)) - corr = cov / np.outer(v, v) + with np.errstate(divide='ignore', invalid='ignore'): + corr = cov / np.outer(v, v) corr[cov == 0] = 0 return corr From d4b0814fdf3e9d6594aed3f1d64d6af2ede28810 Mon Sep 17 00:00:00 2001 From: Richard West Date: Mon, 4 May 2026 09:25:42 -0400 Subject: [PATCH 3/9] Pass scalar temperatures in ESS partition tests The ESS tests passed one-element NumPy arrays into partition-function methods whose Cython signatures expect scalar temperatures. NumPy 1.25 warns about implicit array-to-scalar conversion, and future NumPy releases will make that conversion an error. Defining the test temperature as a float matches the API contract of the methods under test and removes the deprecation warning without changing the expected thermodynamic values. Co-authored-by: Codex with GPT-5.5 --- test/arkane/ess/arkaneGaussianTest.py | 25 ++++++++++++------------- test/arkane/ess/molproTest.py | 11 ++++------- 2 files changed, 16 insertions(+), 20 deletions(-) diff --git a/test/arkane/ess/arkaneGaussianTest.py b/test/arkane/ess/arkaneGaussianTest.py index 8109b708a96..6c2a0bdc632 100644 --- a/test/arkane/ess/arkaneGaussianTest.py +++ b/test/arkane/ess/arkaneGaussianTest.py @@ -93,10 +93,10 @@ def test_load_ethylene_from_gaussian_log_cbsqb3(self): trans = [mode for mode in conformer.modes if isinstance(mode, IdealGasTranslation)][0] rot = [mode for mode in conformer.modes if isinstance(mode, NonlinearRotor)][0] vib = [mode for mode in conformer.modes if isinstance(mode, HarmonicOscillator)][0] - t_list = np.array([298.15], float) - assert abs(trans.get_partition_function(t_list) - 5.83338e6) < 1e1 - assert abs(rot.get_partition_function(t_list) - 2.59622e3) < 1e-2 - assert abs(vib.get_partition_function(t_list) - 1.0481e0) < 1e-4 + t = 298.15 + assert abs(trans.get_partition_function(t) - 5.83338e6) < 1e1 + assert abs(rot.get_partition_function(t) - 2.59622e3) < 1e-2 + assert abs(vib.get_partition_function(t) - 1.0481e0) < 1e-4 assert round(abs(e0 / constants.Na / constants.E_h - -78.467452), 4) == 0 assert conformer.spin_multiplicity == 1 @@ -152,10 +152,10 @@ def test_load_oxygen_from_gaussian_log(self): trans = [mode for mode in conformer.modes if isinstance(mode, IdealGasTranslation)][0] rot = [mode for mode in conformer.modes if isinstance(mode, LinearRotor)][0] vib = [mode for mode in conformer.modes if isinstance(mode, HarmonicOscillator)][0] - t_list = np.array([298.15], float) - assert abs(trans.get_partition_function(t_list) - 7.11169e6) < 1e1 - assert abs(rot.get_partition_function(t_list) - 7.13316e1) < 1e-4 - assert abs(vib.get_partition_function(t_list) - 1.00037e0) < 1e-4 + t = 298.15 + assert abs(trans.get_partition_function(t) - 7.11169e6) < 1e1 + assert abs(rot.get_partition_function(t) - 7.13316e1) < 1e-4 + assert abs(vib.get_partition_function(t) - 1.00037e0) < 1e-4 assert round(abs(e0 / constants.Na / constants.E_h - -150.3784877), 4) == 0 assert conformer.spin_multiplicity == 3 @@ -180,11 +180,10 @@ def test_load_ethylene_from_gaussian_log_g3(self): trans = [mode for mode in conformer.modes if isinstance(mode, IdealGasTranslation)][0] rot = [mode for mode in conformer.modes if isinstance(mode, NonlinearRotor)][0] vib = [mode for mode in conformer.modes if isinstance(mode, HarmonicOscillator)][0] - t_list = np.array([298.15], float) - - assert abs(trans.get_partition_function(t_list) - 5.83338e6) < 1e1 - assert abs(rot.get_partition_function(t_list) - 2.53410e3) < 1e-2 - assert abs(vib.get_partition_function(t_list) - 1.0304e0) < 1e-4 + t = 298.15 + assert abs(trans.get_partition_function(t) - 5.83338e6) < 1e1 + assert abs(rot.get_partition_function(t) - 2.53410e3) < 1e-2 + assert abs(vib.get_partition_function(t) - 1.0304e0) < 1e-4 assert round(abs(e0 / constants.Na / constants.E_h - -78.562189), 4) == 0 assert conformer.spin_multiplicity == 1 diff --git a/test/arkane/ess/molproTest.py b/test/arkane/ess/molproTest.py index 36a121f6238..afeb864fce2 100644 --- a/test/arkane/ess/molproTest.py +++ b/test/arkane/ess/molproTest.py @@ -34,8 +34,6 @@ import os -import numpy as np - import rmgpy.constants as constants from rmgpy.statmech import ( IdealGasTranslation, @@ -129,11 +127,10 @@ def test_load_hosi_from_molpro_log(self): trans = [mode for mode in conformer.modes if isinstance(mode, IdealGasTranslation)][0] rot = [mode for mode in conformer.modes if isinstance(mode, NonlinearRotor)][0] vib = [mode for mode in conformer.modes if isinstance(mode, HarmonicOscillator)][0] - t_list = np.array([298.15], float) - - assert abs(trans.get_partition_function(t_list) - 9.175364e7) < 1e1 - assert abs(rot.get_partition_function(t_list) - 1.00005557e5) < 1e-2 - assert abs(vib.get_partition_function(t_list) - 1.9734989e0) < 1e-4 + t = 298.15 + assert abs(trans.get_partition_function(t) - 9.175364e7) < 1e1 + assert abs(rot.get_partition_function(t) - 1.00005557e5) < 1e-2 + assert abs(vib.get_partition_function(t) - 1.9734989e0) < 1e-4 assert round(abs(e0 / constants.Na / constants.E_h - -768.275662), 4) == 0 assert conformer.spin_multiplicity == 1 From bc743383bf0bdaea321e2fa76ffbd9b398bc7630 Mon Sep 17 00:00:00 2001 From: Richard West Date: Mon, 4 May 2026 09:25:55 -0400 Subject: [PATCH 4/9] Use vector RHS for solvation least squares The solvation K-factor parameter solve used a column-vector right-hand side, causing np.linalg.lstsq to return a two-dimensional parameter array. Converting individual parameters from one-element arrays to floats triggers NumPy scalar-conversion deprecation warnings. Using a one-dimensional right-hand side matches the scalar parameter model and lets NumPy return a one-dimensional solution directly. This removes the warning at the source while leaving the least-squares system and fitted values unchanged. Co-authored-by: Codex with GPT-5.5 --- rmgpy/data/solvation.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/rmgpy/data/solvation.py b/rmgpy/data/solvation.py index bfa1fcb9661..585cbd2e284 100644 --- a/rmgpy/data/solvation.py +++ b/rmgpy/data/solvation.py @@ -2381,7 +2381,7 @@ def get_Kfactor_parameters(self, delG298, delH298, delS298, solvent_name, T_tran # Generate Amatrix and bvector for Ax = b Amatrix = np.zeros((4, 4)) - bvec = np.zeros((4, 1)) + bvec = np.zeros(4) # 1. Tr*ln(K-factor) value at T = 298 K rho_g_298 = get_gas_saturation_density(solvent_name, 298) rho_l_298 = get_liquid_saturation_density(solvent_name, 298) From ce30db6f850d480a6ee0d5da7845a82214b62a78 Mon Sep 17 00:00:00 2001 From: Richard West Date: Mon, 4 May 2026 09:26:06 -0400 Subject: [PATCH 5/9] Adapt Forst optimizer objective to scalar input SciPy optimizers pass the Nelder-Mead search coordinate as an array, but the Forst density-of-states objective is a scalar function. Passing SciPy's array directly into the Cython scalar function caused NumPy array-to-scalar deprecation warnings during pressure-dependence density-of-states calculations. The optimizer now calls a small Python adapter that extracts the scalar value before evaluating phi(), and the returned one-element optimization result is converted explicitly with item(). This keeps the existing optimization method and numerical pathway while making the scalar boundary explicit. Co-authored-by: Codex with GPT-5.5 --- rmgpy/statmech/conformer.pyx | 14 ++++++++++++-- 1 file changed, 12 insertions(+), 2 deletions(-) diff --git a/rmgpy/statmech/conformer.pyx b/rmgpy/statmech/conformer.pyx index a391484fc35..e3e911248db 100644 --- a/rmgpy/statmech/conformer.pyx +++ b/rmgpy/statmech/conformer.pyx @@ -541,6 +541,15 @@ cpdef double phi(double beta, int k, double E, logQ) except -10000000: T = 1.0 / (constants.R * beta) return logQ(T) - k * log(beta) + beta * E + +def _phi_for_optimizer(beta, int k, double E, logQ): + """ + Evaluate :func:`phi` from SciPy optimizers, which pass the scalar search + coordinate as a one-element array. + """ + return phi(float(np.asarray(beta).item()), k, E, logQ) + + @cython.boundscheck(False) @cython.wraparound(False) def get_density_of_states_forst(np.ndarray[np.float64_t, ndim=1] e_list, logQ, int order=1): @@ -557,6 +566,7 @@ def get_density_of_states_forst(np.ndarray[np.float64_t, ndim=1] e_list, logQ, i cdef np.ndarray[np.float64_t, ndim=1] dens_states, sum_states cdef double x, dx, v, E, dE + cdef object opt_result cdef int i, k if order != 1 and order != 2: @@ -580,10 +590,10 @@ def get_density_of_states_forst(np.ndarray[np.float64_t, ndim=1] e_list, logQ, i # Find minimum of phi func x0 arg xtol ftol maxi maxf fullout disp retall callback try: - x = scipy.optimize.fmin(phi, x, (k, E, logQ), 1e-8, 1e-8, 100, 1000, False, False, False, None) + opt_result = scipy.optimize.fmin(_phi_for_optimizer, np.array([x]), (k, E, logQ), 1e-8, 1e-8, 100, 1000, False, False, False, None) except ValueError: break - x = float(x) + x = float(opt_result.item()) dx = 1e-2 * x # Evaluate derivatives needed for steepest descents approximation numerically From c4b34428cf8626858a008bd03a3b34dcd28a45e5 Mon Sep 17 00:00:00 2001 From: Richard West Date: Mon, 4 May 2026 09:26:13 -0400 Subject: [PATCH 6/9] Suppress expected jackknife covariance warning The kinetics family jackknife refit intentionally fits small leave-one-out datasets when estimating rule uncertainty. In those cases SciPy may be unable to estimate parameter covariance, producing an OptimizeWarning even though the fitted rate coefficient is still the value used by the jackknife calculation. The warning suppression is scoped to the jackknife refit that can legitimately produce this warning. This keeps test output focused on actionable warnings without hiding unrelated optimization warnings elsewhere in kinetics fitting. Consideration: this does not improve covariance estimation for sparse jackknife samples; it only acknowledges that covariance is not consumed by this calculation. Co-authored-by: Codex with GPT-5.5 --- rmgpy/data/kinetics/family.py | 22 +++++++++++++++------- 1 file changed, 15 insertions(+), 7 deletions(-) diff --git a/rmgpy/data/kinetics/family.py b/rmgpy/data/kinetics/family.py index 0b43a7134cd..685d92cd778 100644 --- a/rmgpy/data/kinetics/family.py +++ b/rmgpy/data/kinetics/family.py @@ -43,6 +43,7 @@ from copy import deepcopy import numpy as np +from scipy.optimize import OptimizeWarning from sklearn.model_selection import KFold from rmgpy import settings @@ -4654,13 +4655,20 @@ def _make_rule(rr): kin.uncertainty = RateUncertainty(mu=0.0, var=(np.log(fmax) / 2.0) ** 2, N=1, Tref=Tref, data_mean=data_mean, correlation=label) else: if isinstance(rs[0].kinetics, Arrhenius): - dlnks = np.array([ - np.log( - arr().fit_to_reactions(np.delete(rs, i), recipe=recipe) - .to_arrhenius(rxn.get_enthalpy_of_reaction(Tref)) - .get_rate_coefficient(T=Tref) / rxn.get_rate_coefficient(T=Tref) - ) for i, rxn in enumerate(rs) - ]) # 1) fit to set of reactions without the current reaction (k) 2) compute log(kfit/kactual) at Tref + with warnings.catch_warnings(): + # Jackknife refits use sparse samples and only consume fitted rates, not covariance estimates. + warnings.filterwarnings( + "ignore", + message="Covariance of the parameters could not be estimated", + category=OptimizeWarning, + ) + dlnks = np.array([ + np.log( + arr().fit_to_reactions(np.delete(rs, i), recipe=recipe) + .to_arrhenius(rxn.get_enthalpy_of_reaction(Tref)) + .get_rate_coefficient(T=Tref) / rxn.get_rate_coefficient(T=Tref) + ) for i, rxn in enumerate(rs) + ]) # 1) fit to set of reactions without the current reaction (k) 2) compute log(kfit/kactual) at Tref else: dlnks = np.array([ np.log( From d7bde4dff8c7680b39918c425bf1f733862fd3be Mon Sep 17 00:00:00 2001 From: Richard West Date: Mon, 4 May 2026 09:26:19 -0400 Subject: [PATCH 7/9] Guard singular one-state flux solves The pressure-dependence flux filter can encounter a one-isomer steady-state system with no source flux or a zero diagonal matrix entry. Dividing in those singular cases produced runtime warnings and did not yield a meaningful steady-state concentration. Returning None for singular one-state solves reuses the existing fallback path in get_rate_filtered_products(), which already handles unavailable steady-state solutions by falling back to rate-coefficient filtering. This avoids invalid arithmetic while preserving the intended filtering behavior. Consideration: singular one-state networks now take the fallback filtering path instead of returning inf or nan concentrations. Co-authored-by: Codex with GPT-5.5 --- rmgpy/rmg/pdep.py | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/rmgpy/rmg/pdep.py b/rmgpy/rmg/pdep.py index 237cdada2a4..f6daa5292c1 100644 --- a/rmgpy/rmg/pdep.py +++ b/rmgpy/rmg/pdep.py @@ -432,9 +432,14 @@ def solve_ss_network(self, T, P): ind = isomer_spcs.index(self.source[0]) b[ind] = -1.0 # flux at source else: - b = -b / b.sum() # 1.0 flux from source + total_source_flux = b.sum() + if total_source_flux == 0: + return None + b = -b / total_source_flux # 1.0 flux from source if len(b) == 1: + if A[0, 0] == 0: + return None return np.array([b[0] / A[0, 0]]) con = np.linalg.cond(A) From e0db650d6333c455caade14eec9602805cb88185 Mon Sep 17 00:00:00 2001 From: Richard West Date: Mon, 4 May 2026 11:02:28 -0400 Subject: [PATCH 8/9] Render boolean literals in Arkane prettifier The Arkane prettifier is intended to preserve Python literal values while changing layout and whitespace. Rendering True and False as None made the output misleading and caused boolean keyword arguments such as quantum=True and semiclassical=False to lose their meaning in prettified output. Handling bool explicitly in visit_Constant is the right place for this because modern Python ASTs represent booleans as ast.Constant nodes. The check must happen before numeric formatting because bool subclasses int in Python; otherwise True and False could be formatted as 1 and 0. Behavior change: prettify() now renders boolean literals as True and False in keyword arguments and simple list or dict constants. The updated OutputUnitTest expectation records this behavior for the HinderedRotor quantum and semiclassical arguments. Co-authored-by: Codex with GPT-5.5 --- arkane/output.py | 14 +++++++++----- test/arkane/arkaneOutputTest.py | 4 ++-- 2 files changed, 11 insertions(+), 7 deletions(-) diff --git a/arkane/output.py b/arkane/output.py index 1836229c765..1fc8bd1f690 100644 --- a/arkane/output.py +++ b/arkane/output.py @@ -65,8 +65,10 @@ def _is_number(node): return isinstance(node, ast.Constant) and isinstance(node.value, (int, float, complex)) and not isinstance(node.value, bool) @classmethod - def _is_string_or_number(cls, node): - return cls._is_string(node) or cls._is_number(node) + def _is_simple_constant(cls, node): + return cls._is_string(node) or cls._is_number(node) or ( + isinstance(node, ast.Constant) and isinstance(node.value, bool) + ) def visit_Call(self, node): """ @@ -96,7 +98,7 @@ def visit_List(self, node): """ Return a pretty representation of the list represented by `node`. """ - if any([not (self._is_string_or_number(e) or isinstance(e, ast.UnaryOp)) for e in node.elts]): + if any([not (self._is_simple_constant(e) or isinstance(e, ast.UnaryOp)) for e in node.elts]): # Split elements onto multiple lines result = '[\n' self.level += 1 @@ -143,8 +145,8 @@ def visit_Dict(self, node): """ Return a pretty representation of the dict represented by `node`. """ - if (any([not self._is_string_or_number(e) for e in node.keys]) - or any([not self._is_string_or_number(e) for e in node.values])): + if (any([not self._is_simple_constant(e) for e in node.keys]) + or any([not self._is_simple_constant(e) for e in node.values])): # Split elements onto multiple lines result = '{\n' self.level += 1 @@ -167,6 +169,8 @@ def visit_Constant(self, node): """ if isinstance(node.value, str): result = repr(node.value) + elif isinstance(node.value, bool): + result = repr(node.value) elif isinstance(node.value, (int, float, complex)) and not isinstance(node.value, bool): result = '{0:g}'.format(node.value) else: diff --git a/test/arkane/arkaneOutputTest.py b/test/arkane/arkaneOutputTest.py index adde2c6873a..fa1347765d8 100644 --- a/test/arkane/arkaneOutputTest.py +++ b/test/arkane/arkaneOutputTest.py @@ -146,8 +146,8 @@ def test_prettify(self): ], 'kJ/mol', ), - quantum = None, - semiclassical = None, + quantum = True, + semiclassical = False, ), ], spin_multiplicity = 2, From 3ff6418cdebed574b15ca3b4d266c46989b5a227 Mon Sep 17 00:00:00 2001 From: Richard West Date: Wed, 6 May 2026 23:33:51 -0400 Subject: [PATCH 9/9] Return NaN arrays from underdetermined confidence intervals Scalar np.nan is not iterable, so callers that zip or index the returned ci/covariance would raise TypeError when dof <= 0. The efficiency shortcut suggested in code review introduced this error. It looked good in principle, but by committing it through the "accept suggestion" on the PR, we skipped testing if it actually works. It didn't. Anyway, I hope this fixes it. --- arkane/encorr/bac.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/arkane/encorr/bac.py b/arkane/encorr/bac.py index da1cc5d9075..c2e2ccbda07 100644 --- a/arkane/encorr/bac.py +++ b/arkane/encorr/bac.py @@ -1094,7 +1094,7 @@ def get_confidence_intervals(x: np.ndarray, e = y - ypred # Residuals dof = n - p if dof <= 0: - return np.nan, np.nan + return np.full(p, np.nan), np.full((p, p), np.nan) sigma2 = e.T @ weights @ e / dof # MSE cov = sigma2 * np.linalg.inv(x.T @ weights @ x) # covariance matrix se = np.sqrt(np.diag(cov)) # standard error