code: def arima(
      x: np.ndarray,
      order=(0, 0, 0),
      seasonal={"order": (0, 0, 0), "period": 1},
      xreg=None,
      include_mean=True,
      transform_pars=True,
      fixed=None,
      init=None,
      method="CSS",
      SSinit="Gardner1980",
      optim_method="BFGS",
      kappa=1e6,
      tol=1e-8,
      optim_control={"maxiter": 100},
  ):
      SSG = SSinit == "Gardner1980"
      x = x.copy()
  
      def upARIMA(mod, phi, theta):
          p = len(phi)
          q = len(theta)
          mod["phi"] = phi
          mod["theta"] = theta
          r = max(p, q + 1)
          if p > 0:
              mod["T"][:p, 0] = phi
          if r > 1:
              if SSG:
                  mod["Pn"][:r, :r] = getQ0(phi, theta)
              else:
                  raise NotImplementedError('SSinit != "Gardner1980"')
                  # mod['Pn'][:r, :r] = getQ0bis(phi, theta, tol=0)
          else:
              mod["Pn"][0, 0] = 1 / (1 - phi[0] ** 2) if p > 0 else 1
          mod["a"][:] = 0  # a es vector?
          return mod
  
      def arimaSS(y, mod):
          # arima_like(y, phi, theta, delta, a, P, Pn, up, use_resid)
          return arima_like(
              y,
              mod["phi"],
              mod["theta"],
              mod["delta"],
              mod["a"],
              mod["P"],
              mod["Pn"],
              0,
              True,
          )
  
      def armafn(p, x, trans):
          x = x.copy()
          par = coef.copy()
          par[mask] = p
          trarma = arima_transpar(par, arma, trans)
          Z = upARIMA(mod, trarma[0], trarma[1])
          if Z is None:
              return np.finfo(np.float64).max
          if ncxreg > 0:
              x -= np.dot(xreg, par[narma + np.arange(ncxreg)])
          res = arima_like(
              x,
              Z["phi"],
              Z["theta"],
              Z["delta"],
              Z["a"],
              Z["P"],
              Z["Pn"],
              0,
              False,
          )
          if res[2] == 0.0:
              return math.inf
  
          s2 = res[0] / res[2]
          if s2 <= 0:
              return math.nan
          return 0.5 * (math.log(s2) + res[1] / res[2])
  
      def arCheck(ar):
          p = np.argmax(np.append(1, -ar) != 0)
          if not p:
              return True
          coefs = np.append(1, -ar[:p])
          roots = np.polynomial.polynomial.polyroots(coefs)
          return all(np.abs(roots) > 1)
  
      def maInvert(ma):
          q = len(ma)
          q0 = np.argmax(np.append(1, ma) != 0)
          if not q0:
              return ma
          coefs = np.append(1, ma[:q0])
          roots = np.polynomial.polynomial.polyroots(coefs)
          ind = np.abs(roots) < 1
          if any(ind):
              return ma
          if q0 == 1:
              return np.append(1 / ma[0], np.repeat(0, q - q0))
          roots[ind] = 1 / roots[ind]
          x = 1
          for r in roots:
              x = np.append(x, 0) - np.append(0, x) / r
          return x.real[1:], np.repeat(0, q - q0)
  
      if x.ndim > 1:
          raise ValueError(
              "Only implemented for univariate time series"
          )
  
      if x.dtype not in (np.float32, np.float64):
          x = x.astype(np.float64)
      n = len(x)
  
      if len(order) != 3 or any(
          o < 0 or not isinstance(o, int) for o in order
      ):
          raise ValueError(
              f"order must be 3 non-negative integers, got {order}"
          )
      if "order" not in seasonal:
          raise ValueError("order must be a key in seasonal")
      if len(seasonal["order"]) != 3 or any(
          o < 0 or not isinstance(o, int) for o in seasonal["order"]
      ):
          raise ValueError("order must be 3 non-negative integers")
  
      if seasonal["period"] is None or seasonal["period"] == 0:
          warnings.warn("Setting seasonal period to 1")
          seasonal["period"] = 1
  
      # fixed
      # mask
      arma = (
          *order[::2],
          *seasonal["order"][::2],
          seasonal["period"],
          order[1],
          seasonal["order"][1],
      )
      narma = sum(arma[:4])
  
      # xtsp = init x, end x and frequency
      # tsp(x) = None
      Delta = np.array([1.0])
      for i in range(order[1]):
          Delta = tsconv(Delta, np.array([1.0, -1.0]))
  
      for i in range(seasonal["order"][1]):
          Delta = tsconv(
              Delta,
              np.array([1] + [0] * (seasonal["period"] - 1) + [-1]),
          )
      Delta = -Delta[1:]
      nd = order[1] + seasonal["order"][1]
      n_used = (~np.isnan(x)).sum() - len(Delta)
  
      if xreg is None:
          ncxreg = 0
      else:
          if xreg.shape[0] != n:
              raise Exception("lengths of `x` and `xreg` do not match")
  
          if xreg.dtype not in (np.float32, np.float64):
              xreg = xreg.astype(np.float64)
  
          ncxreg = xreg.shape[1]
  
      nmxreg = [f"ex_{i+1}" for i in range(ncxreg)]
      if include_mean and (nd == 0):
          intercept = np.ones(n, dtype=np.float64).reshape(-1, 1)
          if xreg is None:
              xreg = intercept
          else:
              xreg = np.concatenate([intercept, xreg], axis=1)
          ncxreg += 1
          nmxreg = ["intercept"] + nmxreg
  
      # check nas for method CSS-ML
      if method == "CSS-ML":
          anyna = np.isnan(x).any()
          if ncxreg:
              anyna |= np.isnan(xreg).any()
          if anyna:
              method = "ML"
      if method in ["CSS", "CSS-ML"]:
          ncond = order[1] + seasonal["order"][1] * seasonal["period"]
          ncond1 = order[0] + seasonal["order"][0] * seasonal["period"]
          ncond = ncond + ncond1
      else:
          ncond = 0
  
      if fixed is None:
          fixed = np.full(narma + ncxreg, np.nan)
      elif isinstance(fixed, dict):
          add_intercept = "intercept" in nmxreg
          fixed = fixed_params_from_dict(
              fixed,
              order=order,
              seasonal=seasonal,
              intercept=add_intercept,
              n_ex=max(ncxreg - int(add_intercept), 0),
          )
      if len(fixed) != narma + ncxreg:
          raise Exception("wrong length for `fixed`")
      mask = np.isnan(fixed)
  
      no_optim = not mask.any()
  
      if no_optim:
          transform_pars = False
  
      if transform_pars:
          ind = arma[0] + arma[1] + np.arange(arma[2])
          # check masks and more
          if any(~mask[np.arange(arma[0])]) or any(~mask[ind]):
              warnings.warn(
                  "some AR parameters were fixed: setting transform_pars = False"
              )
              transform_pars = False
  
      init0 = np.zeros(narma)
      parscale = np.ones(narma)
  
      # xreg processing
      if ncxreg:
          cn = nmxreg
          orig_xreg = (ncxreg == 1) | (
              ~mask[narma + np.arange(ncxreg)]
          ).any()
          if not orig_xreg:
              _, _, vt = np.linalg.svd(xreg[(~np.isnan(xreg)).all(1)])
              xreg = np.matmul(xreg, vt)
          dx = x
          dxreg = xreg
          if order[1] > 0:
              dx = diff(dx, 1, order[1])
              dxreg = diff(dxreg, 1, order[1])
          if seasonal["period"] > 1 and seasonal["order"][1] > 0:
              dx = diff(dx, seasonal["period"], seasonal["order"][1])
              dxreg = diff(
                  dxreg, seasonal["period"], seasonal["order"][1]
              )
          model = sm.OLS(dx, dxreg)
          result = model.fit()
          fit = {"coefs": result.params, "stderrs": result.bse}
          isna = np.isnan(x) | np.isnan(xreg).any(1)
          n_used = (~isna).sum() - len(Delta)
          init0 = np.append(init0, fit["coefs"])
          ses = fit["stderrs"]
          parscale = np.append(parscale, 10 * ses)
  
      if n_used <= 0:
          raise ValueError("Too few non-missing observations")
  
      if init is not None:
          if len(init) != len(init0):
              raise ValueError(f"init should have length {len(init0)}")
          nan_mask = np.isnan(init)
          if nan_mask.any():
              init[nan_mask] = init0[nan_mask]
          if method == "ML":
              # check stationarity
              if arma[0] > 0:
                  if not arCheck(init[: arma[0]]):
                      raise ValueError("non-stationary AR part")
                  if arma[2] > 0:
                      if not arCheck(
                          init[arma[:2]].sum() + np.arange(arma[2])
                      ):
                          raise ValueError(
                              "non-stationary seasonal AR part"
                          )
                  if transform_pars:
                      init = ARIMA_invtrans(init, arma)
      else:
          init = init0
  
      def arma_css_op(p, x):
          x = x.copy()
          par = coef.copy()
          par[mask] = p
          phi, theta = arima_transpar(par, arma, False)
  
          if ncxreg > 0:
              x -= np.dot(xreg, par[narma + np.arange(ncxreg)])
  
          res, resid = arima_css(x, arma, phi, theta, ncond)
  
          return 0.5 * np.log(res)
  
      coef = np.array(fixed)
      # parscale definition, think about it, scipy doesn't use it
      if method == "CSS":
          if no_optim:
              res = OptimResult(
                  True, 0, np.array([]), 0.0, np.array([])
              )
          else:
              res = minimize(
                  arma_css_op,
                  init[mask],
                  args=(x,),
                  method=optim_method,
                  tol=tol,
                  options=optim_control,
              )
  
          if res.status > 0:
              warnings.warn(
                  f"possible convergence problem: minimize gave code {res.status}]"
              )
          coef[mask] = res.x
          phi, theta = arima_transpar(coef, arma, False)
          mod = make_arima(phi, theta, Delta, kappa)
          if ncxreg > 0:
              x -= np.dot(xreg, coef[narma + np.arange(ncxreg)])
          val = arima_css(x, arma, phi, theta, ncond)
          sigma2 = val[0]
          var = None if no_optim else res.hess_inv / n_used
      else:
          if method == "CSS-ML":
              if not no_optim:
                  res = minimize(
                      arma_css_op,
                      init[mask],
                      args=(x,),
                      method=optim_method,
                      tol=tol,
                      options=optim_control,
                  )
                  # only update the initial parameters if they're valid
                  candidate = init.copy()
                  candidate[mask] = res.x
                  phi, _ = arima_transpar(candidate, arma, False)
                  if np.logical_and(
                      phi > -math.pi / 2, phi < math.pi / 2
                  ).all():
                      init = candidate
              if arma[0] > 0:
                  if not arCheck(init[: arma[0]]):
                      raise ValueError(
                          "non-stationary AR part from CSS"
                      )
              if arma[2] > 0:
                  if not arCheck(
                      init[np.sum(arma[:2])] + np.arange(arma[2])
                  ):
                      raise ValueError(
                          "non-stationary seasonal AR part from CSS"
                      )
              ncond = 0
          if transform_pars:
              init = ARIMA_invtrans(init, arma)
              if arma[1] > 0:
                  ind = arma[0] + np.arange(arma[1])
                  init[ind] = maInvert(init[ind])
              if arma[3] > 0:
                  ind = np.sum(arma[:3]) + np.arange(arma[3])
                  init[ind] = maInvert(init[ind])
          trarma = arima_transpar(init, arma, transform_pars)
          mod = make_arima(trarma[0], trarma[1], Delta, kappa, SSinit)
          if no_optim:
              res = OptimResult(
                  True,
                  0,
                  np.array([]),
                  armafn(np.array([]), x, transform_pars),
                  np.array([]),
              )
          else:
              res = minimize(
                  armafn,
                  init[mask],
                  args=(
                      x,
                      transform_pars,
                  ),
                  method=optim_method,
                  tol=tol,
                  options=optim_control,
              )
          coef[mask] = res.x
          if transform_pars:
              if arma[1] > 0:
                  ind = arma[0] + np.arange(arma[1])
                  if mask[ind].all():
                      coef[ind] = maInvert(coef[ind])
              if arma[3] > 0:
                  ind = np.sum(arma[:3]) + np.arange(arma[3])
                  if mask[ind].all():
                      coef[ind] = maInvert(coef[ind])
              if any(coef[mask] != res.x):
                  oldcode = res.status
                  res = minimize(
                      arma_css_op,
                      coef[mask],
                      args=(x,),
                      method=optim_method,
                      tol=tol,
                      options=optim_control,
                  )
                  res = OptimResult(
                      res.success, oldcode, res.x, res.fun, res.hess_inv
                  )
                  coef[mask] = res.x
              A = arima_gradtrans(coef, arma)
              A = A[np.ix_(mask, mask)]
              sol = np.matmul(res.hess_inv, A) / n_used
              var = np.dot(A, sol)
              coef = arima_undopars(coef, arma)
          else:
              var = None if no_optim else res.hess_inv / n_used
          trarma = arima_transpar(coef, arma, False)
          mod = make_arima(trarma[0], trarma[1], Delta, kappa, SSinit)
          if ncxreg > 0:
              x -= np.dot(xreg, coef[narma + np.arange(ncxreg)])
          val = arimaSS(x, mod)
          val = (val[0], val[3])
          sigma2 = val[0] / n_used
  
      value = 2 * n_used * res.fun + n_used + n_used * np.log(2 * np.pi)
      aic = value + 2 * sum(mask) + 2 if method != "CSS" else np.nan
  
      nm = []
      if arma[0] > 0:
          nm.extend([f"ar{i+1}" for i in range(arma[0])])
      if arma[1] > 0:
          nm.extend([f"ma{i+1}" for i in range(arma[1])])
      if arma[2] > 0:
          nm.extend([f"sar{i+1}" for i in range(arma[2])])
      if arma[3] > 0:
          nm.extend([f"sma{i+1}" for i in range(arma[3])])
      if ncxreg > 0:
          nm += cn
          if not orig_xreg and (var is not None):
              ind = narma + np.arange(ncxreg)
              coef[ind] = np.matmul(vt, coef[ind])
              A = np.identity(narma + ncxreg)
              A[np.ix_(ind, ind)] = vt
              A = A[np.ix_(mask, mask)]
              var = np.matmul(np.matmul(A, var), A.T)
      # if no_optim:
      #     var = pd.DataFrame(var, columns=nm[mask], index=nm[mask])
      resid = val[1]
  
      ans = {
          "coef": dict(zip(nm, coef)),
          "sigma2": sigma2,
          "var_coef": var,
          "mask": mask,
          "loglik": -0.5 * value,
          "aic": aic,
          "arma": arma,
          "residuals": resid,
          #'series': series,
          "code": res.status,
          "n_cond": ncond,
          "nobs": n_used,
          "model": mod,
      }
      return ans
  
imports: import math
  import warnings
  from collections import namedtuple
  from functools import partial
  from typing import Optional, Dict, Union, Tuple
  import numpy
  import pandas
  import statsmodels.api
  from numba import njit
  from scipy.optimize import minimize
  from scipy.stats import norm
  from .mstl import mstl
  from .utils import CACHE, NOGIL
label_desc: Use image node_id 18 for calling a global function with example usage: arima(x, order, seasonal, xreg, include_mean, transform_pars, fixed, init, method, SSinit, optim_method, kappa, tol, optim_control) and returns: ans, mod, arima_like, unknown, all, x,  np, unknown, np, math, math, True, ma, ma, np