<jupyter_start><jupyter_text>Lab03: Logistic Regression.- Student ID: 18127070 - Student name: <NAME> **How to do your homework**You will work directly on this notebook; the word `TODO` indicate the parts you need to do.You can discuss ideas with classmates as well as finding information from the internet, book, etc...; but *this homework must be your*.**How to submit your homework**Before submitting, rerun the notebook (`Kernel` ->` Restart & Run All`).Then create a folder named `ID` (for example, if your ID is 1234567, then name the folder `1234567`). Copy file notebook to this folder, compress and submit it on moodle.**Contents:**- Logistic Regression.  1. Feature Extraction  Import Library<jupyter_code>%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split
from sklearn.metrics import classification_report<jupyter_output><empty_output><jupyter_text>Load data<jupyter_code>from sklearn.datasets import fetch_openml
X, y = fetch_openml('mnist_784', return_X_y=True)
print(X.shape)
print(y.shape)<jupyter_output>(70000, 784)
(70000,)<jupyter_text>Extract Features So we basically have 70000 samples with each sample having 784 features - pixels in this case and a label - the digit the image represent.Let’s play around and see if we can extract any features from the pixels that can be more informative. First I’d like to know more about average intensity - that is the average value of a pixel in an image for the different digits<jupyter_code>labels=np.unique(y)
# print(labels)
n_label=np.unique(y).shape[0]
l_means=np.zeros(shape=n_label,dtype=float) #array stores average intensity for each label

#TODO compute average intensity for each label
cnt = []
for i in range(n_label): 
  cnt.append(i)
for i in range(len(X)):
  idx = int(y[i])
  cnt[idx] += 1 
  l_means[idx] = np.sum(X[i]) / cnt[idx] 
print(l_means)<jupyter_output>[7.94118499 2.18824575 5.54533753 4.8174692  4.28104862 4.2037037
 6.07861087 3.63260274 6.63456754 4.28276159]<jupyter_text>Plot the average intensity using matplotlib<jupyter_code>fig = plt.figure()
ax = fig.add_axes([0,0,1,1])
ax.bar(labels,l_means)<jupyter_output><empty_output><jupyter_text>As we can see there are some differences in intensity. The digit “1” is the less intense while the digit “0” is the most intense. So this new feature seems to have some predictive value if you wanted to know if say your digit is a “1” or not.<jupyter_code>#TODO compute average intensity for each data sample
intensity = np.zeros(shape=y.shape, dtype=float)
for i in range(len(X)):
  intensity[i] = np.mean(X[i])
print(intensity.shape)<jupyter_output>(70000,)<jupyter_text>*Sometimes* people really do not know what are they doing. I am not an exception:<jupyter_code>X_flip=np.flip(X)
symmetry= np.mean((X-X_flip),axis=1)
print(symmetry.shape)<jupyter_output>(70000,)<jupyter_text>I called this feature "symmetry" (though it's not "symmetry" at all). Use visualization method to understand why this thing works. Our new trainning data will have 70000 samples and 2 features: intensity, symmetry.<jupyter_code>#TODO create X_new by horizontal stack intensity and symmetry
X_new = np.column_stack((intensity, symmetry)) #hstack: (140000), vstack(2, 70000)
print(X_new.shape) #it should be (70000,2)<jupyter_output>(70000, 2)<jupyter_text>2. Training  Usually logistic regression is a good first choice for classification. In this homework we use logistic regression for classifying digit 1 images and not digit 1 images.  Normalize dataFirst normalize data using Z-score normalization- **TODO: Study about Z-score normalization**: Z-score normalization is a method of normalizing data that avoids outlier issue. The formula for z-score normalization is: (value - mean) / standard deviation- **TODO: Why should we normalize data?** to avoid outlier issue<jupyter_code>#TODO: normalize X_new
X_new = (X_new - np.mean(X_new)) / np.std(X_new)
print(X_new)<jupyter_output>[[ 0.85342722 -1.6193639 ]
 [ 1.06445354 -0.50547592]
 [ 0.37569199 -1.35218211]
 ...
 [ 0.95427061 -0.19502487]
 [ 0.79632598 -1.04173105]
 [ 1.69918763  0.07215692]]<jupyter_text>Construct  data<jupyter_code>X_new = np.hstack((np.ones((len(X_new), 1)), X_new)) #stack 1s column as usual
y_new=y.astype(int)
y_new[y_new != 1] = 0 # digit 1 -> class 1, other digits -> class 0
y_new=y_new.reshape(-1,1)
print (X_new.shape)
print (y_new.shape)<jupyter_output>(70000, 3)
(70000, 1)<jupyter_text>Split data<jupyter_code>train_X, test_X, train_y, test_y = train_test_split(X_new, y_new, test_size= int(1/3*X.shape[0]))
print(train_X.shape)
print(train_y.shape)<jupyter_output>(46667, 3)
(46667, 1)<jupyter_text>Sigmoid function and derivative of the sigmoid function<jupyter_code>def sigmoid_activation(x):
    """compute the sigmoid activation value for a given input"""
    return 1.0 / (1 + np.exp(-x))
def sigmoid_deriv(x):
    '''compute the derivative of the sigmoid function ASSUMING
    that the input ‘x‘ has already been passed through the sigmoid
    activation function'''
    return x * (1 - x)<jupyter_output><empty_output><jupyter_text>Compute output<jupyter_code>def compute_h(W, X):
    """
    Compute output: Take the dot product between our features ‘X‘ and the weight
    matrix ‘W‘, then pass this value through our sigmoid activation function 
    """
    return sigmoid_activation(X.dot(W))
def predict(W, X):
 
    '''Take the dot product between our features and weight matrix, 
       then pass this value through our sigmoid activation'''
    #........
    preds=sigmoid_activation(X.dot(W))
    # apply a step function to threshold the outputs to binary
    # class labels
    preds[preds <= 0.5] = 0
    preds[preds > 0] = 1

    return preds<jupyter_output><empty_output><jupyter_text>Compute gradient **Loss Function: Average negative log likelihood**$$\mathcal{L}=\dfrac{1}{N} \sum_{i=1}^{N} -\left(y^{i}\ln h_{\mathbf{w}}\left(\mathbf{x}^{i}\right)+\left(1-y^{i}\right)\ln \left(1-h_{\mathbf{w}}\left(x^{i}\right)\right)\right) $$$$\text{Sigmoid Activation: } z= \sigma \left(h\right)= \dfrac{1}{1+e^{-h}}$$$$\text{Cross-entropy: } J(w)=-\left({ylog(z)+(1-y)log(1-z)}\right)$$$$\text{Chain rule: } \dfrac{\partial J(w)}{\partial w}=\dfrac{\partial J(w)}{\partial z} \dfrac{\partial z}{\partial h}\dfrac{\partial h}{\partial w}  $$$$\dfrac{\partial J(w)}{\partial z}=-\left(\dfrac{y}{z}-\dfrac{1-y}{1-z}\right)=\dfrac{z-y}{z(1-z)}$$$$\dfrac{\partial z}{\partial h}=z(1-z)$$$$\dfrac{\partial h}{\partial w}=X$$$$\dfrac{\partial J(w)}{\partial w}=X^T(z-y)$$<jupyter_code>def compute_gradient(h, error):
    """
    This is the gradient descent update of "average negative loglikelihood" loss function. 
    In lab02 our loss function is "sum squared error".
    """
    #TODO
    gradient = train_X.T.dot(error) / len(train_X)
    return gradient
def train(W,train_X, train_y, learning_rate, num_epochs, losses):
    for epoch in np.arange(0, num_epochs):
        h=compute_h(W,train_X)
        error = h - train_y
        loss = np.mean(- train_y * np.log(h) - (1 - train_y) * np.log(1 - h))
        losses.append(loss)
        gradient=compute_gradient(h, error)
        W += -learning_rate * gradient
        if ((epoch+1)%1000==0): print ('Epoch %d, loss %.3f' %(epoch+1, loss))
        
    return W<jupyter_output><empty_output><jupyter_text>Train our model<jupyter_code>W = np.random.randn(train_X.shape[1], 1)
losses=[]
num_epochs=40000
learning_rate=0.01
W=train(W,train_X, train_y, learning_rate, num_epochs , losses)
x_preds=predict(W ,train_X)
train_err = np.mean(x_preds != train_y) * 100
print ('=' * 50)
print ('Train err of final w: ', train_err)<jupyter_output>Epoch 1000, loss 0.295
Epoch 2000, loss 0.259
Epoch 3000, loss 0.242
Epoch 4000, loss 0.233
Epoch 5000, loss 0.228
Epoch 6000, loss 0.224
Epoch 7000, loss 0.221
Epoch 8000, loss 0.219
Epoch 9000, loss 0.217
Epoch 10000, loss 0.215
Epoch 11000, loss 0.214
Epoch 12000, loss 0.213
Epoch 13000, loss 0.212
Epoch 14000, loss 0.211
Epoch 15000, loss 0.211
Epoch 16000, loss 0.210
Epoch 17000, loss 0.209
Epoch 18000, loss 0.209
Epoch 19000, loss 0.208
Epoch 20000, loss 0.208
Epoch 21000, loss 0.207
Epoch 22000, loss 0.207
Epoch 23000, loss 0.207
Epoch 24000, loss 0.206
Epoch 25000, loss 0.206
Epoch 26000, loss 0.206
Epoch 27000, loss 0.205
Epoch 28000, loss 0.205
Epoch 29000, loss 0.205
Epoch 30000, loss 0.204
Epoch 31000, loss 0.204
Epoch 32000, loss 0.204
Epoch 33000, loss 0.204
Epoch 34000, loss 0.203
Epoch 35000, loss 0.203
Epoch 36000, loss 0.203
Epoch 37000, loss 0.203
Epoch 38000, loss 0.203
Epoch 39000, loss 0.203
Epoch 40000, loss 0.202
=================================================[...]<jupyter_text>3. Evaluate our model<jupyter_code>preds = predict(W, train_X)
print(classification_report(train_y, preds))
preds = predict(W, test_X)
print(classification_report(test_y, preds))<jupyter_output>precision    recall  f1-score   support

           0       0.93      0.98      0.95     20670
           1       0.74      0.41      0.53      2663

    accuracy                           0.92     23333
   macro avg       0.84      0.70      0.74     23333
weighted avg       0.91      0.92      0.91     23333<jupyter_text>**TODO: Comment on the result**: Because the percentage of loss of our training is quite low, so the classification table for both train and test are quite similar about result<jupyter_code><jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12754237288135592, 'fraction_numerical': 0.0847457627118644, 'mean_word_length': 4.47304347826087, 'pattern_counts': {'":': 0, '<': 57, '<?xml version=': 0, '>': 60, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '16', 'dir': 'jupyter-structured-clean-dedup', 'id': '88383', 'license': 'MIT', 'n_tokens_mistral': 3880, 'n_tokens_neox': 3267, 'n_words': 998, 'path': 'Lab03_LogisticRegression.ipynb', 'repo_name': 'IceIce1ce/Lab-Logistic-Regression-HCMUS', 'size': '9446'}