docs: def _asym_transform_deriv_shape(systematic_utilities,
                                  alt_IDs,
                                  rows_to_alts,
                                  eta,
                                  ref_position=None,
                                  dh_dc_array=None,
                                  fill_dc_d_eta=None,
                                  output_array=None,
                                  *args, **kwargs):
      """
      Parameters
      ----------
      systematic_utilities : 1D ndarray.
          Contains the systematic utilities for each each available alternative
          for each observation. All elements should be ints, floats, or longs.
      alt_IDs : 1D ndarray.
         All elements should be ints. There should be one row per obervation per
         available alternative for the given observation. Elements denote the
         alternative corresponding to the given row of the design matrix.
      rows_to_alts : 2D ndarray.
          There should be one row per observation per available alternative and
          one column per possible alternative. This matrix maps the rows of the
          design matrix to the possible alternatives for this dataset.
      eta : 1D ndarray.
          Each element should be an int, float, or long. There should be one
          value per transformed shape parameter. Note that if there are J
          possible alternatives in the dataset, then there should be J - 1
          elements in `eta`.
      ref_position : int.
          Specifies the position in the array of natural shape parameters that
          should be equal to 1 - the sum of the other elements. Specifies the
          alternative in the ordered array of unique alternatives that is not
          having its shape parameter estimated (to ensure identifiability).
      dh_dc_array : 2D scipy sparse matrix.
          Its data is to be replaced with the correct derivatives of the
          transformed index vector with respect to the shape parameter vector.
          Should have shape
          `(systematic_utilities.shape[0], rows_to_alts.shape[1])`.
      fill_dc_d_eta : callable.
          Should accept `eta` and `ref_position` and return a 2D numpy array
          containing the derivatives of the 'natural' shape parameter vector with
          respect to the vector of transformed shape parameters.
      output_array : 2D numpy matrix.
          This matrix's data is to be replaced with the correct derivatives of
          the transformed systematic utilities with respect to the vector of
          transformed shape parameters. Should have shape
          `(systematic_utilities.shape[0], shape_params.shape[0])`.
  
      Returns
      -------
      output_array : 2D ndarray.
          The shape of the returned array will be
          `(systematic_utilities.shape[0], shape_params.shape[0])`. The returned
          array specifies the derivative of the transformed utilities with
          respect to the shape parameters.
      """
      ##########
      # Convert the reduced shape parameters to the natural shape parameters
      ##########
      natural_shape_params = _convert_eta_to_c(eta, ref_position)
  
      ##########
      # Calculate the derivative of the transformed utilities with respect to
      # the vector of natural shape parameters, c
      ##########
      # Create a vector which contains the appropriate shape for each row in the
      # design matrix. Note as long as natural_shape_params is a numpy array,
      # then long_shapes will be a numpy array.
      long_shapes = rows_to_alts.dot(natural_shape_params)
  
      # Calculate d_ln(long_shape)_d_long_shape
      d_lnShape_dShape = 1.0 / long_shapes
      # Guard against overflow
      d_lnShape_dShape[np.isposinf(d_lnShape_dShape)] = max_comp_value
  
      # Calculate d_ln((1-long_shape)/(J-1))_d_long_shape
      d_lnShapeComp_dShape = -1.0 / (1 - long_shapes)
      # Guard against overflow
      d_lnShapeComp_dShape[np.isneginf(d_lnShapeComp_dShape)] = -max_comp_value
  
      # Differentiate the multiplier with respect to natural_shape_j.
      deriv_multiplier = ((systematic_utilities >= 0) * d_lnShape_dShape +
                          (systematic_utilities < 0) * d_lnShapeComp_dShape)
      # assert not np.isnan(deriv_multiplier).any()
  
      # Calculate the derivative of h_ij with respect to natural_shape_j.
      # Store these derivatives in their respective places in the dh_dc array
      # Note that d_hij_d_ck = 0 for k != j
      dh_dc_values = d_lnShape_dShape - systematic_utilities * deriv_multiplier
      # Guard against overflow
      dh_dc_values[np.isinf(dh_dc_values)] = -1 * max_comp_value
      # Assign the computed values to the scipy sparse array
      dh_dc_array.data = dh_dc_values
  
      ##########
      # Calculate the derivative of the natural shape parameters, c with
      # respect to the vector of reduced shape parameters, eta
      ##########
      # Return the matrix of dh_d_eta. Note the matrix should be of dimension
      # (systematic_utilities.shape[0], shape_params.shape[0])
      # Note the calculation is essentially dh_dc * dc_d_eta = dh_d_eta
      output_array[:, :] = dh_dc_array.dot(fill_dc_d_eta(natural_shape_params,
                                                         ref_position))
      return output_array