input: Write a python function to implement gradient descent for linear regression and compute the cost at each iteration.
  assert np.allclose(gradient_descent(np.array([[1, 1], [1, 2], [1, 3]]), np.array([1, 2, 3]), np.array([0, 0]), 0.1, 10)[0], [0.26, 0.46], atol=0.01)
  assert np.allclose(gradient_descent(np.array([[1, 1], [1, 2], [1, 3], [1, 4], [1, 5]]), np.array([1, 2, 3, 4, 5]), np.array([0, 0]), 0.01, 1000)[0], [0.00, 1.00], atol=0.01)
  assert np.allclose(gradient_descent(np.array([[1, 1], [1, 2], [1, 3]]), np.array([2, 2.5, 3.5]), np.array([0, 0]), 0.1, 100)[0], [0.5, 1.0], atol=0.01)
output: import numpy as np
  
  def gradient_descent(X, y, theta, alpha, num_iters):
      m = len(y)
      J_history = []
  
      for i in range(num_iters):
          prediction = X.dot(theta)
          error = prediction - y
          gradient = (1/m) * X.T.dot(error)
          theta = theta - alpha * gradient
          cost = (1/(2*m)) * np.sum(error**2)
          J_history.append(cost)
  
      return theta, J_history