Fix type annotations, RuntimeWarning, and half-life display issues

README.md

@@ -100,14 +100,15 @@ Enforced automatically during estimation.
 
 ## API Reference
 
-### `estimate_gjr_garch_x(returns, exog_vars, method='SLSQP', verbose=False)`
+### `estimate_gjr_garch_x(returns, exog_vars, method='SLSQP', max_iter=1000, verbose=False)`
 
 Main estimation function.
 
 **Parameters:**
 - `returns`: `pd.Series` of log returns (recommend × 100 for numerical stability)
 - `exog_vars`: `pd.DataFrame` of exogenous variables, aligned with returns index
 - `method`: Optimization method (`'SLSQP'`, `'L-BFGS-B'`, `'trust-constr'`)
+- `max_iter`: Maximum number of optimizer iterations (default: 1000)
 - `verbose`: Print estimation progress
 
 **Returns:** `GJRGARCHXResults` object

src/gjr_garch_x/__init__.py

@@ -36,7 +36,7 @@
 
 import warnings
 from dataclasses import dataclass
-from typing import Dict, List, Optional, Tuple
+from typing import Any, Dict, List, Optional, Tuple, cast
 
 import numpy as np
 import pandas as pd
@@ -48,6 +48,13 @@
 __author__ = "Murad Farzulla"
 __email__ = "murad@farzulla.org"
 
+# Constants
+# Threshold for detecting near-unit-root persistence in GARCH models.
+# When persistence is >= 0.9999, the half-life calculation becomes numerically
+# unstable and produces misleading negative values. This threshold prevents
+# displaying nonsensical half-life values for processes very close to non-stationary.
+NEAR_UNIT_ROOT_THRESHOLD = 0.9999
+
 __all__ = [
     "estimate_gjr_garch_x",
     "GJRGARCHXResults",
@@ -114,7 +121,7 @@ class GJRGARCHXResults:
     iterations: int
     n_obs: int = 0
 
-    def __post_init__(self):
+    def __post_init__(self) -> None:
         """Set n_obs from volatility length if not provided."""
         if self.n_obs == 0 and len(self.volatility) > 0:
             object.__setattr__(self, "n_obs", len(self.volatility))
@@ -166,9 +173,11 @@ def summary(self) -> str:
         lines.append("")
         lines.append(f"Persistence (α + β + |γ|/2): {persistence:.4f}")
 
-        if 0 < persistence < 1:
+        if 0 < persistence < NEAR_UNIT_ROOT_THRESHOLD:
             half_life = -np.log(0.5) / np.log(persistence)
             lines.append(f"Half-life of shocks:         {half_life:.1f} periods")
+        elif persistence >= NEAR_UNIT_ROOT_THRESHOLD:
+            lines.append("Half-life of shocks:         ∞ (near unit root)")
 
         # Unconditional variance (if stationary)
         omega = self.params.get("omega", 0)
@@ -281,14 +290,14 @@ def __init__(
     def _unpack_params(self, params: np.ndarray) -> Dict[str, float]:
         """Unpack parameter vector into named dictionary."""
         param_dict = {
-            "omega": params[0],
-            "alpha": params[1],
-            "gamma": params[2],
-            "beta": params[3],
-            "nu": params[4],
+            "omega": float(params[0]),
+            "alpha": float(params[1]),
+            "gamma": float(params[2]),
+            "beta": float(params[3]),
+            "nu": float(params[4]),
         }
         for i, name in enumerate(self.exog_names):
-            param_dict[name] = params[5 + i]
+            param_dict[name] = float(params[5 + i])
         return param_dict
 
     def _variance_recursion(
@@ -312,7 +321,7 @@ def _variance_recursion(
 
         variance = np.zeros(self.n_obs)
         mean_return = self.returns.mean()
-        residuals = (self.returns - mean_return).values
+        residuals = np.asarray((self.returns - mean_return).values)
 
         # Initialize with unconditional variance estimate
         variance[0] = np.var(self.returns)
@@ -368,7 +377,7 @@ def _log_likelihood(self, params: np.ndarray) -> float:
         except (ValueError, OverflowError, RuntimeWarning):
             return 1e8
 
-    def _parameter_constraints(self) -> List[Dict]:
+    def _parameter_constraints(self) -> List[Dict[str, Any]]:
         """Define optimization constraints including stationarity."""
         return [
             {"type": "ineq", "fun": lambda x: x[0] - 1e-8},  # omega > 0
@@ -394,8 +403,10 @@ def _get_starting_values(self) -> np.ndarray:
             5.0,  # nu
         ])
         if self.has_exog:
-            start_vals = np.append(start_vals, np.zeros(self.n_exog))
-        return start_vals
+            exog_zeros = np.zeros(self.n_exog)
+            result = np.concatenate([start_vals, exog_zeros])
+            return cast(np.ndarray, result)
+        return cast(np.ndarray, start_vals)
 
     def estimate(
         self,
@@ -425,12 +436,12 @@ def estimate(
 
         start_vals = self._get_starting_values()
 
-        bounds = [
+        bounds: List[Tuple[Optional[float], Optional[float]]] = [
             (1e-8, None),  # omega > 0
             (1e-8, 0.3),  # alpha
             (-0.5, 0.5),  # gamma (leverage)
             (1e-8, 0.95),  # beta
-            (2.1, 50),  # nu
+            (2.1, 50.0),  # nu
         ]
         # Exogenous coefficients are unbounded
         for _ in range(self.n_exog):
@@ -533,14 +544,19 @@ def _compute_standard_errors(
         try:
             hessian = self._numerical_hessian(params)
             cov_matrix = np.linalg.inv(hessian)
-            std_errs = np.sqrt(np.diag(cov_matrix))
+            diag = np.diag(cov_matrix)
+
+            # Check for negative variances before sqrt (ill-conditioned Hessian)
+            if np.any(diag < 0):
+                nan_dict = {name: np.nan for name in self.param_names}
+                return (nan_dict, nan_dict)
+
+            std_errs = np.sqrt(diag)
 
             dof = self.n_obs - self.n_params
             if dof <= 0:
-                return (
-                    {name: np.nan for name in self.param_names},
-                    {name: np.nan for name in self.param_names},
-                )
+                nan_dict = {name: np.nan for name in self.param_names}
+                return (nan_dict, nan_dict)
 
             t_stats = params / std_errs
             pvals = 2 * (1 - student_t.cdf(np.abs(t_stats), dof))
@@ -551,15 +567,13 @@ def _compute_standard_errors(
             )
 
         except (np.linalg.LinAlgError, ValueError):
-            return (
-                {name: np.nan for name in self.param_names},
-                {name: np.nan for name in self.param_names},
-            )
+            nan_dict = {name: np.nan for name in self.param_names}
+            return (nan_dict, nan_dict)
 
     def _numerical_hessian(self, params: np.ndarray, h: float = 1e-5) -> np.ndarray:
         """Compute numerical Hessian via central differences."""
         n = len(params)
-        hessian = np.zeros((n, n))
+        hessian: np.ndarray = np.zeros((n, n))
 
         for i in range(n):
             for j in range(n):