In other words, the more “linear” data is, the more precise we can predict.

Like this:

In multidimensional space: $y = a_1 \cdot x_1 + a_2\cdot x_2 + ….. + b$ (depend on dimension, this can be line - 2D or plane - 3D or hyperplane - xD).

Here, b is also called bias value. We all know in linear algebra, the line is perfectly start from the origin of the 2D coordinate O(0, 0), but in reality, there is always a bias value cause the line a little bit far from the origin.

In the picture above, y-intercept (or bias) is 5.

Therefore, we obtain loss function (loss here means “error”, because the distance from point to fitting line is error, we want to minimize the error).

• N is number of data points (number of blue dot in the above graph).
• Inside the sum is the Error formula. Because we assume the perfect line is starting from the origin, which have function: y’ = w . X and the reality line is y = ax + b —> we just simply implement: error = y - y’, it may be negative, we must use square or absolute value. However, we have to take the derivative afterward, so just use square (or can say “Euclidean Norm).
• X bar is a vector with multiple columns, in multidimensional data, X bar is a matrix with corresponding columns and rows.

• 1/2 is a scalar for easy to compute, because in reality, the dataset is massive. When taking the derivative, the 0.5 scalar is automatically reduced.

  Another way to represent loss function:

  \[L(w) = \frac{1}{2} \cdot \|y - \bar{X} \cdot w \|_2^2\]

However, in many problems, the function is hard to find the derivative or hard to compute the value for equation: derivative = 0. Luckily, in linear regression, it is possible.

By using partial derivative, obtain:

• For the derivative = 0, obtain:

• So, the equation have solutions if the vector / matrix multiplication in the lefthand side is invertible (non singular). We obtain the solution for w:

  \[w = (\bar{X}^T \cdot \bar{X})^{-1} \cdot \bar{X}^T \cdot y\]

  The complex thing: $(\bar{X}^T\cdot\bar{X})^{-1}\cdot \bar{X}^T$ is called “Pseudo inverse” (just the name). And the equation above is called “Normal Equation”.

• Work well with small number of features, when the features scale bigger, the normal equation is slow (features are different from data point, features mean the category of data, for example, for house price prediction, features are area, number of bedrooms, quality of neighborhood and so on).
• The X Transpose * X must be invertible, otherwise, must use regularization to handle
• Time complexity is O(n^3) where n is number of features.

To be more specific:

Firstly, we can see here, it is easy for us to do this as a programmer 😉, obviously, a for loop might help!

Secondly, in detail, the complex part is just the derivative of the loss function, which we consider here is linear regression.

Thirdly, the “eta” symbol is called “learning rate”. Clearly, machine is not human, it can not walk on it own, we need to set the pace of its walk (to the valley), the pace also reflects the speed of this algorithms. If learning rate is too big (like moving a big step downward), the algorithm will never give precise answer (because, once again, it is a machine, not human). And vice versa, too small —> too slow to reach the valley (or even never reach).

—> Overall, imagine a man called “Loss”, he carried a bunch of package on his body, he want to go down to the valley as soon as possible, so he take off all the necessary things to obtain a lightweight body (called “gradient”), we set for him an initial place, which is x0, and require him to go with, say, 0.1 kilometer at each step (each loop).

—> In conclusion, learning rate eta is too low —> slow, too high —> never reach optimal point (valley). Remark:

• In the formula, we update x1 by let x0 - eta * derivative of loss function and so on for other values —> do it in one loop with identified iteration (we set iteration).
• For different x0, we often get nearly the same solution

Example:

$f(x) = x^2 + 5 \cdot sin(x)$

$\Rightarrow f’(x) = 2 \cdot x + 5 \cdot cos(x)$

  import numpy as np
  import matplotlib.pyplot as plt
        
  def Grad(x):
    return 2*x+ 5*np.cos(x)
        
  def Cost(x):
    return x**2 + 5*np.sin(x)
          
  def GD(x_0, weight):
    res = []
    for i in range(100):
      x_0 = x_0 - weight*Grad(x_0)
      res.append(x_0)
      if abs(Grad(x_0)) < 1e-3:
        break
    return (res[-1], i)

Testing:

  print(GD(-5, 0.1))
  print(GD(5, 0.1))

—> the minimum value obtained after 11 loops (x_0 = -5). The same with 29 loops (x_0 = 5)

Plot:

  x = np.arange(-5, 5, 1) #x_initial = -5
  f = Cost(x)
  point_x = GD(-5, 0.1)[0]
  point_y = Cost(point_x)
  print(mark)
  plt.plot(x, f)
  plt.scatter(point_x, point_y)
  plt.xlabel("x")
  plt.ylabel("f(x)")
  plt.show()

For linear regression:

  def Loss(weight):
  	return 0.5/N * np.linalg.norm((y - np.dot(X, weight))) ** 2
  def Gradient(weight):
  	(1 / N) * np.dot(X.T, ( np.dot(X, weight) - y))
  def Update(weight_init, eta, iteration):
  	train_size = len(X)
  	weight = weight_init
  	for i in range(iteration):
  		loss = Loss(weight)
  		weight = weight - (eta / train_size) * gradient(weight)
  		loss /= train_size
  	return weight, loss

—> weights are the output of the function.

• Easy to compute (use the power of machine), high adaptability.
• Matrix of dataset doesn’t have to be invertible, square,…
• Complexity: O(n^2)

So, the formula for 1 point is:

To be easier, just using 2 loops, each epoch (outside loop) runs n times (n is number of datapoints), means that each epoch, our algorithm will take a look to a whole data based on 1 datapoint (given in random order). During an epoch, shuffle data of X parameter is require.

  def stochastic_gradient_descent(X, y, w, eta, iter):
      # batch size = 1
      train_loss = []
      temp = np.random.randint(N)
      for epoch in range(iter):
          np.random.shuffle(X)
          for i in range(N):
              xi = X[temp : temp + 1]
              yi = y[temp : temp + 1]
              l = (np.dot(xi, w) - yi)**2
              grad = 2* xi.T.dot(np.dot(xi, w) - yi)
              w = w - (eta) * grad
              train_loss.append(l[0][0]/N)
      return w, train_loss

• faster in the speed of convergence —> can find the optimal point quickly.
• require a small number of iteration.

Now the formula is

• n is number of batch size = number of points we want to take. A good tip is just taking a number that is a perfect square root i.e 2^2, 2^3, 2^4,….
• x is the data from i to i+n.
• Like SGD, shuffle X is required to have a better results 😂

  def minibatch_gradient_descent(X, y, w, eta, iter, batch_size):
      train_loss = []
      for epoch in range(iter):
          np.random.shuffle(X)
          for i in range(0, features, batch_size):
              X_i = X[i:i+batch_size-1]
              y_i = y[i:i+batch_size-1]
              l = (np.dot(X_i, w) - y_i)**2
              grad = 2* X_i.T.dot(np.dot(X_i, w) - y_i)  
              w = w - (eta) * grad
              train_loss.append(l[-1][-1]/N)
      return w, train_loss

• Faster in speed
• Small number of iterations
• Useful with a huge number of data.

Vietnamese:

https://machinelearningcoban.com/2017/01/16/gradientdescent2/

English:

https://cs231n.github.io/neural-networks-3/

https://ruder.io/optimizing-gradient-descent/index.html#stochasticgradientdescent

https://en.wikipedia.org/wiki/Stochastic_gradient_descent