code: def _incremental_mean_and_var(
      X, last_mean, last_variance, last_sample_count
  ):
      """Calculate mean update and a Youngs and Cramer variance update.
  
      last_mean and last_variance are statistics computed at the last step by the
      function. Both must be initialized to 0.0. In case no scaling is required
      last_variance can be None. The mean is always required and returned because
      necessary for the calculation of the variance. last_n_samples_seen is the
      number of samples encountered until now.
  
      From the paper "Algorithms for computing the sample variance: analysis and
      recommendations", by Chan, Golub, and LeVeque.
  
      Parameters
      ----------
      X : array-like, shape (n_samples, n_features)
          Data to use for variance update
  
      last_mean : array-like, shape: (n_features,)
  
      last_variance : array-like, shape: (n_features,)
  
      last_sample_count : array-like, shape (n_features,)
  
      Returns
      -------
      updated_mean : array, shape (n_features,)
  
      updated_variance : array, shape (n_features,)
          If None, only mean is computed
  
      updated_sample_count : array, shape (n_features,)
  
      Notes
      -----
      NaNs are ignored during the algorithm.
  
      References
      ----------
      T. Chan, G. Golub, R. LeVeque. Algorithms for computing the sample
          variance: recommendations, The American Statistician, Vol. 37, No. 3,
          pp. 242-247
  
      Also, see the sparse implementation of this in
      `utils.sparsefuncs.incr_mean_variance_axis` and
      `utils.sparsefuncs_fast.incr_mean_variance_axis0`
      """
      # old = stats until now
      # new = the current increment
      # updated = the aggregated stats
      last_sum = last_mean * last_sample_count
      new_sum = _safe_accumulator_op(np.nansum, X, axis=0)
  
      new_sample_count = np.sum(~np.isnan(X), axis=0)
      updated_sample_count = last_sample_count + new_sample_count
  
      updated_mean = (last_sum + new_sum) / updated_sample_count
  
      if last_variance is None:
          updated_variance = None
      else:
          new_unnormalized_variance = (
              _safe_accumulator_op(np.nanvar, X, axis=0)
              * new_sample_count
          )
          last_unnormalized_variance = last_variance * last_sample_count
  
          with cupyx.errstate(divide=None, invalid=None):
              last_over_new_count = last_sample_count / new_sample_count
              updated_unnormalized_variance = (
                  last_unnormalized_variance
                  + new_unnormalized_variance
                  + last_over_new_count
                  / updated_sample_count
                  * (last_sum / last_over_new_count - new_sum) ** 2
              )
  
          zeros = last_sample_count == 0
          updated_unnormalized_variance[
              zeros
          ] = new_unnormalized_variance[zeros]
          updated_variance = (
              updated_unnormalized_variance / updated_sample_count
          )
  
      return updated_mean, updated_variance, updated_sample_count
  
imports: from cuml.internals.safe_imports import gpu_only_import_from
  from cuml.internals.safe_imports import gpu_only_import
label_desc: Use image node_id 2 for calling a global function with example usage: _incremental_mean_and_var(X, last_mean, last_variance, last_sample_count) and returns: updated_mean,  updated_variance,  updated_sample_count