docs: def minimize(self, theta_init, max_iter=50, callback=None, disp=0, tau=(10., 2., 2.), tol=1e-3):
          """
          Minimize a list of objectives using a proximal consensus algorithm
  
          Parameters
          ----------
          theta_init : ndarray
              Initial parameter vector (numpy array)
  
          max_iter : int, optional
              Maximum number of iterations to run (default: 50)
  
          callback : function, optional
              a function that gets called on each iteration with the following arguments: the current parameter
              value (ndarray), and a dictionary that contains a information about the status of the algorithm
  
          disp : int, optional
              determines how much information to display when running. Ranges from 0 (nothing) to 3 (lots of information)
  
          Returns
          -------
          theta : ndarray
              The parameters found after running the optimization procedure
  
          Other Parameters
          ----------------
          tau : (float, float, float), optional
              initial, increment and decrement parameters for the momentum scheduler (default: (10, 2, 2))
  
          tol : float, optional
              residual tolerance for assessing convergence. if both the primal and dual residuals are less
              than this value, then the algorithm has converged (default: 1e-3)
          """
  
          # get list of objectives for this parameter
          num_obj = len(self.objectives)
          assert num_obj >= 1, "There must be at least one objective!"
  
          # initialize lists of primal and dual variable copies, one for each objective
          orig_shape = theta_init.shape
          primals = [theta_init.flatten() for _ in range(num_obj)]
          duals = [np.zeros(theta_init.size) for _ in range(num_obj)]
          theta_avg = np.mean(primals, axis=0).ravel()
  
          # initialize penalty parameter
          tau = namedtuple('tau', ('init', 'inc', 'dec'))(*tau)
          rho = tau.init
  
          # store cumulative runtimes of each iteration, starting now
          tstart = time.time()
  
          # clear metadata
          self.metadata = defaultdict(list)
  
          # run ADMM iterations
          self.converged = False
          for cur_iter in range(max_iter):
  
              # store the parameters from the previous iteration
              theta_prev = theta_avg
  
              # update each primal variable copy by taking a proximal step via each objective
              for varidx, dual in enumerate(duals):
                  primals[varidx] = self.objectives[varidx]((theta_prev - dual).reshape(orig_shape), rho).ravel()
  
              # average primal copies
              theta_avg = np.mean(primals, axis=0)
  
              # update the dual variables (after primal update has finished)
              for varidx, primal in enumerate(primals):
                  duals[varidx] += primal - theta_avg
  
              # compute primal and dual residuals
              primal_resid = float(np.sum([np.linalg.norm(primal - theta_avg) for primal in primals]))
              dual_resid = num_obj * rho ** 2 * np.linalg.norm(theta_avg - theta_prev)
  
              # update penalty parameter according to primal and dual residuals
              # (see sect. 3.4.1 of the Boyd and Parikh ADMM paper)
              if primal_resid > tau.init * dual_resid:
                  rho *= float(tau.inc)
              elif dual_resid > tau.init * primal_resid:
                  rho /= float(tau.dec)
  
              # update metadata for this iteration
              self.metadata['Primal resid'].append(primal_resid)
              self.metadata['Dual resid'].append(dual_resid)
              self.metadata['Time (s)'].append(time.time() - tstart)
              self.metadata['rho'].append(rho)
  
              # invoke the callback function with the current parameters and
              # history
              if callback is not None:
  
                  # get the metadata from this iteration
                  data = valmap(last, self.metadata)
  
                  callback(theta_avg.reshape(orig_shape), data)
  
              # update the display
              self.update_display(cur_iter + 1, disp)
  
              # check for convergence
              if (primal_resid <= tol) & (dual_resid <= tol):
                  self.converged = True
                  break
  
          # clean up display
          self.update_display(-1, disp)
  
          # store and return final parameters
          self.theta = theta_avg.reshape(orig_shape)
          return self.theta