<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def sample(self, X, y, quantity='y', sample_at_X=None, weights=None, n_draws=100, n_bootstraps=5, objective='auto'):
 """Simulate from the posterior of the coefficients and smoothing params. Samples are drawn from the posterior of the coefficients and smoothing parameters given the response in an approximate way. The GAM must already be fitted before calling this method; if the model has not been fitted, then an exception is raised. Moreover, it is recommended that the model and its hyperparameters be chosen with `gridsearch` (with the parameter `keep_best=True`) before calling `sample`, so that the result of that gridsearch can be used to generate useful response data and so that the model's coefficients (and their covariance matrix) can be used as the first bootstrap sample. These samples are drawn as follows. Details are in the reference below. 1. ``n_bootstraps`` many "bootstrap samples" of the response (``y``) are simulated by drawing random samples from the model's distribution evaluated at the expected values (``mu``) for each sample in ``X``. 2. A copy of the model is fitted to each of those bootstrap samples of the response. The result is an approximation of the distribution over the smoothing parameter ``lam`` given the response data ``y``. 3. Samples of the coefficients are simulated from a multivariate normal using the bootstrap samples of the coefficients and their covariance matrices. Notes ----- A ``gridsearch`` is done ``n_bootstraps`` many times, so keep ``n_bootstraps`` small. Make ``n_bootstraps < n_draws`` to take advantage of the expensive bootstrap samples of the smoothing parameters. Parameters X : array of shape (n_samples, m_features) empirical input data y : array of shape (n_samples,) empirical response vector quantity : {'y', 'coef', 'mu'}, default: 'y' What quantity to return pseudorandom samples of. If `sample_at_X` is not None and `quantity` is either `'y'` or `'mu'`, then samples are drawn at the values of `X` specified in `sample_at_X`. sample_at_X : array of shape (n_samples_to_simulate, m_features) or None, optional Input data at which to draw new samples. Only applies for `quantity` equal to `'y'` or to `'mu`'. If `None`, then `sample_at_X` is replaced by `X`. weights : np.array of shape (n_samples,) sample weights n_draws : positive int, optional (default=100) The number of samples to draw from the posterior distribution of the coefficients and smoothing parameters n_bootstraps : positive int, optional (default=5) The number of bootstrap samples to draw from simulations of the response (from the already fitted model) to estimate the distribution of the smoothing parameters given the response data. If `n_bootstraps` is 1, then only the already fitted model's smoothing parameter is used, and the distribution over the smoothing parameters is not estimated using bootstrap sampling. objective : string, optional (default='auto' metric to optimize in grid search. must be in ['AIC', 'AICc', 'GCV', 'UBRE', 'auto'] if 'auto', then grid search will optimize GCV for models with unknown scale and UBRE for models with known scale. Returns ------- draws : 2D array of length n_draws Simulations of the given `quantity` using samples from the posterior distribution of the coefficients and smoothing parameter given the response data. Each row is a pseudorandom sample. If `quantity == 'coef'`, then the number of columns of `draws` is the number of coefficients (`len(self.coef_)`). Otherwise, the number of columns of `draws` is the number of rows of `sample_at_X` if `sample_at_X` is not `None` or else the number of rows of `X`. References Simon N. Wood, 2006. Generalized Additive Models: an introduction with R. Section 4.9.3 (pages 198–199) and Section 5.4.2 (page 256–257). """

        if quantity not in {'mu', 'coef', 'y'}:
            raise ValueError("`quantity` must be one of 'mu', 'coef', 'y';"
                             " got {}".format(quantity))

        coef_draws = self._sample_coef(
            X, y, weights=weights, n_draws=n_draws,
            n_bootstraps=n_bootstraps, objective=objective)

        if quantity == 'coef':
            return coef_draws

        if sample_at_X is None:
            sample_at_X = X

        linear_predictor = self._modelmat(sample_at_X).dot(coef_draws.T)
        mu_shape_n_draws_by_n_samples = self.link.mu(
            linear_predictor, self.distribution).T
        if quantity == 'mu':
            return mu_shape_n_draws_by_n_samples
        else:
            return self.distribution.sample(mu_shape_n_draws_by_n_samples)