text: import numpy as np
  from prml.linear.regression import Regression
  
  
  class BayesianRegression(Regression):
      """
      Bayesian regression model
  
      w ~ N(w|0, alpha^(-1)I)
      y = X @ w
      t ~ N(t|X @ w, beta^(-1))
      """
  
      def __init__(self, alpha:float=1., beta:float=1.):
          self.alpha = alpha # w의 prior 분산과 관련
          self.beta = beta # 회귀모형의 epsilon의 분산과 괸련
          self.w_mean = None # w의 posterior mean이 바인딩될 곳
          self.w_precision = None # w의 posterior variance가 바인딩될 곳
  
      def _is_prior_defined(self) -> bool:
          return self.w_mean is not None and self.w_precision is not None
  
      def _get_prior(self, ndim:int) -> tuple:
          if self._is_prior_defined():
              return self.w_mean, self.w_precision
          else:
              return np.zeros(ndim), self.alpha * np.eye(ndim) # 객체에서 처음 fit을 호출하면 초기의 prior를 준다
  
      def fit(self, X:np.ndarray, t:np.ndarray):
          """
          bayesian update of parameters given training dataset
  
          Parameters
          ----------
          X : (N, n_features) np.ndarray
              training data independent variable
          t : (N,) np.ndarray
              training data dependent variable
          """
  
          mean_prev, precision_prev = self._get_prior(np.size(X, 1)) # 이전의 posterior는 현재의 prior
  
          # posterior 구하기
          w_precision = precision_prev + self.beta * X.T @ X 
          w_mean = np.linalg.solve(
              w_precision,
              precision_prev @ mean_prev + self.beta * X.T @ t
          )
          self.w_mean = w_mean
          self.w_precision = w_precision
          self.w_cov = np.linalg.inv(self.w_precision)
  
      def predict(self, X:np.ndarray, return_std:bool=False, sample_size:int=None):
          """
          return mean (and standard deviation) of predictive distribution
  
          Parameters
          ----------
          X : (N, n_features) np.ndarray
              independent variable
          return_std : bool, optional
              flag to return standard deviation (the default is False)
          sample_size : int, optional
              number of samples to draw from the predictive distribution
              (the default is None, no sampling from the distribution)
  
          Returns
          -------
          y : (N,) np.ndarray
              mean of the predictive distribution
          y_std : (N,) np.ndarray
              standard deviation of the predictive distribution
          y_sample : (N, sample_size) np.ndarray
              samples from the predictive distribution
          """
  
          if sample_size is not None: # 1. sampling mode
              w_sample = np.random.multivariate_normal(
                  self.w_mean, self.w_cov, size=sample_size
              )
              y_sample = X @ w_sample.T
              return y_sample
          y = X @ self.w_mean # predictive mean
          if return_std: # 2. distribution mode (mean & variance)
              y_var = 1 / self.beta + np.sum(X @ self.w_cov * X, axis=1) # predictive variance
              y_std = np.sqrt(y_var)
              return y, y_std
          return y
  
  # TODO January 04, 2020: 챕터 3.3을 공부하고 주석달기