code: def posterior(self, x, sigma=1.):
          
          pr0 = 1. / self.scale**2  # prior precision
          prd = x.size / sigma**2  # data precision
          varp = 1. / (pr0 + prd)  # posterior variance
          mu = varp * (pr0 * self.loc + prd * x.mean())
          return Normal(loc=mu, scale=np.sqrt(varp))
docs: Model is X_1,...,X_n ~ N(theta, sigma^2), theta~self, sigma fixed