code: def _calc_deriv_c_with_respect_to_eta(natural_shapes,
                                        ref_position,
                                        output_array=None):
      
      # Generate a list of the indices which indicate the columns to be
      # selected from a 2D numpy array of
      # np.diag(natural_shapes) - np.outer(natural_shapes, natural_shapes)
      columns_to_be_kept = range(natural_shapes.shape[0])
      columns_to_be_kept.remove(ref_position)
  
      # Calculate and store the derivative of the natural shape parameters
      # with respect to the reduced shape parameters.
      output_array[:, :] = (np.diag(natural_shapes) -
                            np.outer(natural_shapes,
                                     natural_shapes))[:, columns_to_be_kept]
  
      return output_array
docs: Parameters
      ----------
      natural_shapes : 1D ndarray.
          Should have one element per available alternative in the dataset whose
          choice situations are being modeled. Should have at least
          `ref_position` elements in it.
      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 (in order to ensure
          identifiability).
      output_array : 2D ndarray.
          This array is to have its data overwritten with the correct derivatives
          of the natural shape parameters with respect to transformed shape
          parameters. Should have shape ==
          `(natural_shapes.shape[0], natural_shapes.shape[0] - 1)`.
  
      Returns
      -------
      output_array : 2D ndarray.
          Has shape == (natural_shapes.shape[0], natural_shapes.shape[0] - 1).
          Will contain the derivative of the shape parameters, with
          respect to the underlying 'transformed' shape parameters.