<jupyter_start><jupyter_text>4. Training Models Treating algorithms as black boxes is definitely faster and it does the trick, but understand what is under the hood is essential to really have a firm grasp on our systems and being able to adapt them to what we need.  In this section, we will :1. Work on a Linear Regression model, training it using both a **closed-form** equation and **Gradient Descent**2. Introduce **Polynomial regression**, useful for non-linear datasets but more prone to overfitting (regularization will be needed)3. We will look at two more models commonly used for classification: **Logistic Regression** and **Softmax Regression**  Linear regressionGenerally, a linear model makes a prediction by simply computing a weighted sum of the **input features**, plus a constant called the **bias term** (also called the _intercept term_).$y = \theta_0 + \theta_1x_1 + \theta_2x_2 + \cdots + \theta_nx_n$ Or, in vectorized form: $y = h_{\theta}(x) = \theta \cdot x$ To train this model, what we want to do is **minimize the Mean Square Error (MSE)**: $\text{MSE }(X, h_\theta) = \frac{1}{m} \sum^{m}_{i = 1}(\theta^T x^{(i)} - y^{(i)})^2$  Normal EquationIt is simply the _closed-form solution_ (an equation that computes the result directly without intermediate steps) for our minimization problem:  $ \theta_{min} = (X^T X)^{-1} X^T y$  Example:<jupyter_code>import numpy as np

X = 2 * np.random.rand(100, 1)
y = 4 + 3 * X + np.random.randn(100, 1)
X_b = np.c_[np.ones((100, 1)), X] # add x0 = 1 to each instance
theta_best = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)
theta_best
# making predictions using theta_best

X_new = np.array([[0], [2]])
X_new_b = np.c_[np.ones((2, 1)), X_new] # add x0 = 1 to each instance
y_predict = X_new_b.dot(theta_best)
y_predict
from matplotlib import pyplot as plt

plt.plot(X_new, y_predict, "r-")
plt.plot(X, y, "b.")
plt.axis([0, 2, 0, 15])
plt.show()<jupyter_output><empty_output><jupyter_text>Now let's use Scikit-learn:<jupyter_code>from sklearn.linear_model import LinearRegression

lin_reg = LinearRegression()
lin_reg.fit(X, y)
lin_reg.intercept_, lin_reg.coef_
lin_reg.predict(X_new)<jupyter_output><empty_output><jupyter_text>Gradient Descent Gradient Descent is a very generic optimization algorithm capable of iteratively finding optimal solutions to a wide range of problems. Given its iterative approach, GD may run into issues when it comes to local/global minima. Luckily for us, the Linear Regression cost function is _convex_ (it only has one global minimum).    **Note**: when using Gradient Descent, we should ensure that all features have a similar scale for faster convergence.   Batch Gradient Descent It's the most common implementation, going over the whole training set at each GD step.  $ \theta_{\text{(next step)}} = \theta - \eta \nabla_{\theta} \text{MSE}(\theta)$ Python implementation:<jupyter_code>eta = 0.1 # learning rate
n_iterations = 1000
m = 100

theta = np.random.randn(2,1) # random initialization

for iteration in range(n_iterations):
    gradients = 2/m * X_b.T.dot(X_b.dot(theta) - y)
    theta = theta - eta * gradients
theta<jupyter_output><empty_output><jupyter_text>Same as our Normal Equation result above!  Stochastic Gradient DescentObviously, Batch Gradient Descent can get very slow for big datasets, since it is working on the whole training set. On the opposite of the spectrum, we have Stochastic GD, using a random instance in the training set at every step. Due to its stochastic nature, it has a better chance to find global optima for functions which are not strictly convex, but it is also not bound to settle for any specific minimum, but rather staying in a "good enough" area.A solution to this dilemma may be to gradually decrease the learning rate (according to a **learning schedule**).<jupyter_code>n_epochs = 50
t0, t1 = 5, 50 # learning schedule hyperparameters

def learning_schedule(t):
    return t0 / (t + t1)

theta = np.random.randn(2,1) # random initialization

for epoch in range(n_epochs):
    for i in range(m):
        random_index = np.random.randint(m)
        xi = X_b[random_index:random_index + 1]
        yi = y[random_index:random_index + 1]
        gradients = 2 * xi.T.dot(xi.dot(theta) - yi)
        eta = learning_schedule(epoch * m + i)
        theta = theta - eta * gradients
theta<jupyter_output><empty_output><jupyter_text>Pretty good, but not optimal. **Note**: When using Stochastic Gradient Descent, the training instances must be independent and identically distributed (IID), to ensure that the parameters get pulled towards the global optimum, on average. A simple way to ensure this is to **shuffle the instances during training**. Example:<jupyter_code>from sklearn.linear_model import SGDRegressor # Stochastic Gradient Descent Regressor

sgd_reg = SGDRegressor(max_iter=1000, tol=1e-3, penalty=None, eta0=0.1)
sgd_reg.fit(X, y.ravel())
sgd_reg.intercept_, sgd_reg.coef_<jupyter_output><empty_output><jupyter_text>Mini-batch Gradient Descent The middle way. Instead of going full-training-set (Batch) or one-instance-only (Stochastic) Mini-batch GD computes the gradients on small random sets of instances called _... drum rolls..._ mini-batches.To conclude, here is a summary table for the different algorithms available for Linear Regression: ![Comparison Table of algorithms for Linear Regression](images/4.Lin_Reg_algos_comparison.jpg)  Polynomial RegressionWe can also use linear models to fit nonlinear data by adding powers of each feature as new features, then train a linear model on this extended set of features. This technique is called **Polynomial Regression**. Example:<jupyter_code># generating non-linear (quadratic) data

m = 100
X = 6 * np.random.rand(m, 1) - 3
y = 0.5 * X**2 + X + 2 + np.random.randn(m, 1)<jupyter_output><empty_output><jupyter_text>We can now use Scikit-Learn’s `PolynomialFeatures` class to transform our training data, adding the square of each feature in the training set as new features.<jupyter_code>from sklearn.preprocessing import PolynomialFeatures

poly_features = PolynomialFeatures(degree=2, include_bias=False) # excluding bias term
X_poly = poly_features.fit_transform(X)
X[0]
X_poly[0] # original feature + squared feature<jupyter_output><empty_output><jupyter_text>We can now fit a Linear Regression to this extended training data:<jupyter_code>lin_reg = LinearRegression()
lin_reg.fit(X_poly, y)
lin_reg.intercept_, lin_reg.coef_<jupyter_output><empty_output><jupyter_text>**Note**: Given **n** features, `PolynomialFeatures(degree=d)` returns an array of $\displaystyle \frac{(n + d)!}{d!n!}$  Learning CurvesA standard way to understand if our model is under/overfitting the data is to look at the learning curves: plots of the model’s performance on the training set and the validation set as a function of the training set size (or the training iteration). Generally speaking:1. Converging + high error = **underfitting**2. Not converging + lower error = **overfitting** (adding more data can help) **Note**: A model’s generalization error can be expressed as the sum of three different errors:1. **Bias**: due to wrong assumptions > likely to underfitting2. **Variance**: due to excessive sensitivity to small variations in the training data > likely to overfit3. **Irreducible error**: due to the intrinsic noisiness of the data > we can reduce it by cleaning up the data (broken sources, outliers, etc.)  Regularized Linear ModelsFor a linear model, regularization is typically achieved by constraining the weights of the model.Let's look at three different ways on how this can be achieved：  Ridge RegressionA regularization term $\alpha \sum_{i = 1}^n \theta^2_i$ is added to the cost function we are trying to minimize: $J(\theta) = MSE(\theta) + \alpha\frac{1}{2}\sum^n_{i=1} \theta^2_i$ **Note 1**: the regularization term should only be added to the cost function during training.  **Note 2**: Remember to scale the data before performing regularization.  Lasso RegressionShorter for _Least Absolute Shrinkage and Selection Operator Regression_. It uses the $l_1$ norm of the weight vector instead of half the square of the $l_2$ norm: $J(\theta) = MSE(\theta) + \alpha\sum^n_{i=1} | \theta_i |$ An important characteristic of Lasso Regression is that it tends to completely eliminate the weights of the least important features (i.e., set them to zero).  In other words, Lasso Regression automatically performs feature selection and outputs a **sparse model** (i.e., with few nonzero feature weights).  Elastic NetIt's a mix of Ridge and Lasso, which we can control with a parameter $r$: * r = 0 > Ridge Regression* r = 1 > Lasso Regression$J(\theta) = MSE(\theta) + r\alpha\sum^n_{i=1} |\theta_i| + \frac{1-r}{2} \alpha \sum^n_{i=1} \theta^2_i$  How to choose a Linear ModelTwo general tips:1. **Ridge** is a good default, but if you suspect that only a few features are actually useful, you should prefer Lasso or Elastic Net2. **Elastic Net** is preferred over Lasso since Lasso may behave erratically when the number of features is greater than the number of training instances or when several features are strongly correlated  Early stoppingWe can also decide to stop GD as soon as it reaches a certain error level. This may be particularly helpful in complex models where the **validation error** stops decreasing and actually starts to go back up.  For Stochastic or Mini-Batch GD, since they are not _smooth_, we may want to stop only after the error has been under the minimimum for a certain number of epochs.   Logistic Regression Logistic Regression (also called _Logit Regression_) is commonly used to estimate the probability that an instance belongs to a particular class.  If the probability is above a certain **threshold**, the model predicts that the instance belongs to that class making Logistic Regression a **binary classifier**.  A Logistic Regression model computes a weighted sum of the input features (plus a bias term), but instead of outputting the result directly like the Linear Regression model does, it outputs the **logistic** of this result: $p = h_{\theta}(x) = \sigma(x^T \theta)$Where $\sigma$ is a sigmoid function which outputs a number between 0 and 1:$\displaystyle \sigma(t) = \frac{1}{1+ exp(-t)}$  Cost Function Here is our _log loss_:$J(\theta) = -\frac{1}{m} \sum^m_{i=1} [y^{(i)} log(p^{(i)}) + (1-y^{(i)}) log (1-p^{(i)})]$  Decision BoundariesLet's build a classifier to detect the Iris-Virginica using a well-known dataset and Logistic Regression:<jupyter_code>import numpy as np
from sklearn import datasets
import matplotlib.pyplot as plt

iris = datasets.load_iris()
list(iris.keys())
X = iris["data"][:, 3:] # petal width
y = (iris["target"] == 2).astype(np.int) # 1 if Iris-Virginica, else 0
# training our Logistic Regression
from sklearn.linear_model import LogisticRegression

log_reg = LogisticRegression(solver='lbfgs')
log_reg.fit(X, y)
# estimated probs

X_new = np.linspace(0, 3, 1000).reshape(-1, 1)
y_proba = log_reg.predict_proba(X_new)
plt.plot(X_new, y_proba[:, 1], "g-", label="Iris-Virginica")
plt.plot(X_new, y_proba[:, 0], "b--", label="Not Iris-Virginica")<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.09596500624888413, 'fraction_numerical': 0.012319228709159078, 'mean_word_length': 5.0198817839871035, 'pattern_counts': {'":': 0, '<': 37, '<?xml version=': 0, '>': 42, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '9', 'dir': 'jupyter-structured-clean-dedup', 'id': '245373', 'license': 'MIT', 'n_tokens_mistral': 3469, 'n_tokens_neox': 3280, 'n_words': 1408, 'path': '.ipynb_checkpoints/04.Training_Models-checkpoint.ipynb', 'repo_name': 'GiuseppeBaldini/ML-Scikit-Keras-TensorFlow', 'size': '11202'}