adding material to week 42
This commit is contained in:
Morten Hjorth-Jensen
File diff suppressed because one or more lines are too long
@@ -0,0 +1,441 @@
|
||||
TITLE: Exercises Week 42: Logistic Regression and Optimization, reminders from week 38 and week 40
|
||||
AUTHOR: Morten Hjorth-Jensen {copyright, 1999-present|CC BY-NC} at Department of Physics and Center for Computing in Science Education, University of Oslo & Department of Physics and Astronomy and Facility for Rare Isotope Beams, Michigan State University
|
||||
DATE: October 14-18, 2024
|
||||
|
||||
!split
|
||||
===== The logistic function =====
|
||||
|
||||
A widely studied model, is the
|
||||
perceptron model, which is an example of a ``hard classification'' model. We
|
||||
have used this model when we discussed neural networks as
|
||||
well. Each datapoint is deterministically assigned to a category (i.e
|
||||
$y_i=0$ or $y_i=1$). In many cases it is favorable to have a ``soft''
|
||||
classifier that outputs the probability of a given category rather
|
||||
than a single value. For example, given $x_i$, the classifier
|
||||
outputs the probability of being in a category $k$. Logistic regression
|
||||
is the most common example of a so-called soft classifier. In logistic
|
||||
regression, the probability that a data point $x_i$
|
||||
belongs to a category $y_i=\{0,1\}$ is given by the so-called logit function (or Sigmoid) which is meant to represent the likelihood for a given event,
|
||||
!bt
|
||||
\[
|
||||
p(t) = \frac{1}{1+\mathrm \exp{-t}}=\frac{\exp{t}}{1+\mathrm \exp{t}}.
|
||||
\]
|
||||
!et
|
||||
Note that $1-p(t)= p(-t)$.
|
||||
|
||||
!split
|
||||
===== Examples of likelihood functions used in logistic regression and nueral networks =====
|
||||
|
||||
|
||||
The following code plots the logistic function, the step function and other functions we will encounter from here and on.
|
||||
|
||||
|
||||
!bc pycod
|
||||
"""The sigmoid function (or the logistic curve) is a
|
||||
function that takes any real number, z, and outputs a number (0,1).
|
||||
It is useful in neural networks for assigning weights on a relative scale.
|
||||
The value z is the weighted sum of parameters involved in the learning algorithm."""
|
||||
|
||||
import numpy
|
||||
import matplotlib.pyplot as plt
|
||||
import math as mt
|
||||
|
||||
z = numpy.arange(-5, 5, .1)
|
||||
sigma_fn = numpy.vectorize(lambda z: 1/(1+numpy.exp(-z)))
|
||||
sigma = sigma_fn(z)
|
||||
|
||||
fig = plt.figure()
|
||||
ax = fig.add_subplot(111)
|
||||
ax.plot(z, sigma)
|
||||
ax.set_ylim([-0.1, 1.1])
|
||||
ax.set_xlim([-5,5])
|
||||
ax.grid(True)
|
||||
ax.set_xlabel('z')
|
||||
ax.set_title('sigmoid function')
|
||||
|
||||
plt.show()
|
||||
|
||||
"""Step Function"""
|
||||
z = numpy.arange(-5, 5, .02)
|
||||
step_fn = numpy.vectorize(lambda z: 1.0 if z >= 0.0 else 0.0)
|
||||
step = step_fn(z)
|
||||
|
||||
fig = plt.figure()
|
||||
ax = fig.add_subplot(111)
|
||||
ax.plot(z, step)
|
||||
ax.set_ylim([-0.5, 1.5])
|
||||
ax.set_xlim([-5,5])
|
||||
ax.grid(True)
|
||||
ax.set_xlabel('z')
|
||||
ax.set_title('step function')
|
||||
|
||||
plt.show()
|
||||
|
||||
"""tanh Function"""
|
||||
z = numpy.arange(-2*mt.pi, 2*mt.pi, 0.1)
|
||||
t = numpy.tanh(z)
|
||||
|
||||
fig = plt.figure()
|
||||
ax = fig.add_subplot(111)
|
||||
ax.plot(z, t)
|
||||
ax.set_ylim([-1.0, 1.0])
|
||||
ax.set_xlim([-2*mt.pi,2*mt.pi])
|
||||
ax.grid(True)
|
||||
ax.set_xlabel('z')
|
||||
ax.set_title('tanh function')
|
||||
|
||||
plt.show()
|
||||
!ec
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
!split
|
||||
===== Two parameters =====
|
||||
|
||||
We assume now that we have two classes with $y_i$ either $0$ or $1$. Furthermore we assume also that we have only two parameters $\beta$ in our fitting of the Sigmoid function, that is we define probabilities
|
||||
!bt
|
||||
\begin{align*}
|
||||
p(y_i=1|x_i,\bm{\beta}) &= \frac{\exp{(\beta_0+\beta_1x_i)}}{1+\exp{(\beta_0+\beta_1x_i)}},\nonumber\\
|
||||
p(y_i=0|x_i,\bm{\beta}) &= 1 - p(y_i=1|x_i,\bm{\beta}),
|
||||
\end{align*}
|
||||
!et
|
||||
where $\bm{\beta}$ are the weights we wish to extract from data, in our case $\beta_0$ and $\beta_1$.
|
||||
|
||||
Note that we used
|
||||
!bt
|
||||
\[
|
||||
p(y_i=0\vert x_i, \bm{\beta}) = 1-p(y_i=1\vert x_i, \bm{\beta}).
|
||||
\]
|
||||
!et
|
||||
|
||||
|
||||
!split
|
||||
===== The cost function =====
|
||||
|
||||
Reordering the logarithms, we can rewrite the _cost/loss_ function as
|
||||
!bt
|
||||
\[
|
||||
\mathcal{C}(\bm{\beta}) = \sum_{i=1}^n \left(y_i(\beta_0+\beta_1x_i) -\log{(1+\exp{(\beta_0+\beta_1x_i)})}\right).
|
||||
\]
|
||||
!et
|
||||
|
||||
The maximum likelihood estimator is defined as the set of parameters that maximize the log-likelihood where we maximize with respect to $\beta$.
|
||||
Since the cost (error) function is just the negative log-likelihood, for logistic regression we have that
|
||||
!bt
|
||||
\[
|
||||
\mathcal{C}(\bm{\beta})=-\sum_{i=1}^n \left(y_i(\beta_0+\beta_1x_i) -\log{(1+\exp{(\beta_0+\beta_1x_i)})}\right).
|
||||
\]
|
||||
!et
|
||||
This equation is known in statistics as the _cross entropy_. Finally, we note that just as in linear regression,
|
||||
in practice we often supplement the cross-entropy with additional regularization terms, usually $L_1$ and $L_2$ regularization as we did for Ridge and Lasso regression.
|
||||
|
||||
!split
|
||||
===== Minimizing the cross entropy =====
|
||||
|
||||
The cross entropy is a convex function of the weights $\bm{\beta}$ and,
|
||||
therefore, any local minimizer is a global minimizer.
|
||||
|
||||
|
||||
Minimizing this
|
||||
cost function with respect to the two parameters $\beta_0$ and $\beta_1$ we obtain
|
||||
|
||||
!bt
|
||||
\[
|
||||
\frac{\partial \mathcal{C}(\bm{\beta})}{\partial \beta_0} = -\sum_{i=1}^n \left(y_i -\frac{\exp{(\beta_0+\beta_1x_i)}}{1+\exp{(\beta_0+\beta_1x_i)}}\right),
|
||||
\]
|
||||
!et
|
||||
and
|
||||
!bt
|
||||
\[
|
||||
\frac{\partial \mathcal{C}(\bm{\beta})}{\partial \beta_1} = -\sum_{i=1}^n \left(y_ix_i -x_i\frac{\exp{(\beta_0+\beta_1x_i)}}{1+\exp{(\beta_0+\beta_1x_i)}}\right).
|
||||
\]
|
||||
!et
|
||||
|
||||
!split
|
||||
===== A more compact expression =====
|
||||
|
||||
Let us now define a vector $\bm{y}$ with $n$ elements $y_i$, an
|
||||
$n\times p$ matrix $\bm{X}$ which contains the $x_i$ values and a
|
||||
vector $\bm{p}$ of fitted probabilities $p(y_i\vert x_i,\bm{\beta})$. We can rewrite in a more compact form the first
|
||||
derivative of cost function as
|
||||
|
||||
!bt
|
||||
\[
|
||||
\frac{\partial \mathcal{C}(\bm{\beta})}{\partial \bm{\beta}} = -\bm{X}^T\left(\bm{y}-\bm{p}\right).
|
||||
\]
|
||||
!et
|
||||
|
||||
If we in addition define a diagonal matrix $\bm{W}$ with elements
|
||||
$p(y_i\vert x_i,\bm{\beta})(1-p(y_i\vert x_i,\bm{\beta})$, we can obtain a compact expression of the second derivative as
|
||||
|
||||
!bt
|
||||
\[
|
||||
\frac{\partial^2 \mathcal{C}(\bm{\beta})}{\partial \bm{\beta}\partial \bm{\beta}^T} = \bm{X}^T\bm{W}\bm{X}.
|
||||
\]
|
||||
!et
|
||||
|
||||
!split
|
||||
===== Extending to more predictors =====
|
||||
|
||||
Within a binary classification problem, we can easily expand our model to include multiple predictors. Our ratio between likelihoods is then with $p$ predictors
|
||||
!bt
|
||||
\[
|
||||
\log{ \frac{p(\bm{\beta}\bm{x})}{1-p(\bm{\beta}\bm{x})}} = \beta_0+\beta_1x_1+\beta_2x_2+\dots+\beta_px_p.
|
||||
\]
|
||||
!et
|
||||
Here we defined $\bm{x}=[1,x_1,x_2,\dots,x_p]$ and $\bm{\beta}=[\beta_0, \beta_1, \dots, \beta_p]$ leading to
|
||||
!bt
|
||||
\[
|
||||
p(\bm{\beta}\bm{x})=\frac{ \exp{(\beta_0+\beta_1x_1+\beta_2x_2+\dots+\beta_px_p)}}{1+\exp{(\beta_0+\beta_1x_1+\beta_2x_2+\dots+\beta_px_p)}}.
|
||||
\]
|
||||
!et
|
||||
|
||||
!split
|
||||
===== Including more classes =====
|
||||
|
||||
Till now we have mainly focused on two classes, the so-called binary
|
||||
system. Suppose we wish to extend to $K$ classes. Let us for the sake
|
||||
of simplicity assume we have only two predictors. We have then following model
|
||||
|
||||
!bt
|
||||
\[
|
||||
\log{\frac{p(C=1\vert x)}{p(K\vert x)}} = \beta_{10}+\beta_{11}x_1,
|
||||
\]
|
||||
!et
|
||||
and
|
||||
!bt
|
||||
\[
|
||||
\log{\frac{p(C=2\vert x)}{p(K\vert x)}} = \beta_{20}+\beta_{21}x_1,
|
||||
\]
|
||||
!et
|
||||
and so on till the class $C=K-1$ class
|
||||
!bt
|
||||
\[
|
||||
\log{\frac{p(C=K-1\vert x)}{p(K\vert x)}} = \beta_{(K-1)0}+\beta_{(K-1)1}x_1,
|
||||
\]
|
||||
!et
|
||||
|
||||
and the model is specified in term of $K-1$ so-called log-odds or
|
||||
_logit_ transformations.
|
||||
|
||||
|
||||
!split
|
||||
===== More classes =====
|
||||
|
||||
In our discussion of neural networks we will encounter the above again
|
||||
in terms of a slightly modified function, the so-called _Softmax_ function.
|
||||
|
||||
The softmax function is used in various multiclass classification
|
||||
methods, such as multinomial logistic regression (also known as
|
||||
softmax regression), multiclass linear discriminant analysis, naive
|
||||
Bayes classifiers, and artificial neural networks. Specifically, in
|
||||
multinomial logistic regression and linear discriminant analysis, the
|
||||
input to the function is the result of $K$ distinct linear functions,
|
||||
and the predicted probability for the $k$-th class given a sample
|
||||
vector $\bm{x}$ and a weighting vector $\bm{\beta}$ is (with two
|
||||
predictors):
|
||||
|
||||
!bt
|
||||
\[
|
||||
p(C=k\vert \mathbf {x} )=\frac{\exp{(\beta_{k0}+\beta_{k1}x_1)}}{1+\sum_{l=1}^{K-1}\exp{(\beta_{l0}+\beta_{l1}x_1)}}.
|
||||
\]
|
||||
!et
|
||||
It is easy to extend to more predictors. The final class is
|
||||
!bt
|
||||
\[
|
||||
p(C=K\vert \mathbf {x} )=\frac{1}{1+\sum_{l=1}^{K-1}\exp{(\beta_{l0}+\beta_{l1}x_1)}},
|
||||
\]
|
||||
!et
|
||||
|
||||
and they sum to one.
|
||||
|
||||
|
||||
!split
|
||||
===== Wisconsin Cancer Data =====
|
||||
|
||||
We show here how we can use a simple regression case on the breast
|
||||
cancer data using Logistic regression as our algorithm for
|
||||
classification.
|
||||
|
||||
|
||||
!bc pycod
|
||||
import matplotlib.pyplot as plt
|
||||
import numpy as np
|
||||
from sklearn.model_selection import train_test_split
|
||||
from sklearn.datasets import load_breast_cancer
|
||||
from sklearn.linear_model import LogisticRegression
|
||||
|
||||
# Load the data
|
||||
cancer = load_breast_cancer()
|
||||
|
||||
X_train, X_test, y_train, y_test = train_test_split(cancer.data,cancer.target,random_state=0)
|
||||
print(X_train.shape)
|
||||
print(X_test.shape)
|
||||
# Logistic Regression
|
||||
logreg = LogisticRegression(solver='lbfgs')
|
||||
logreg.fit(X_train, y_train)
|
||||
print("Test set accuracy with Logistic Regression: {:.2f}".format(logreg.score(X_test,y_test)))
|
||||
!ec
|
||||
|
||||
!split
|
||||
===== Using the correlation matrix =====
|
||||
|
||||
In addition to the above scores, we could also study the covariance (and the correlation matrix).
|
||||
We use _Pandas_ to compute the correlation matrix.
|
||||
!bc pycod
|
||||
import matplotlib.pyplot as plt
|
||||
import numpy as np
|
||||
from sklearn.model_selection import train_test_split
|
||||
from sklearn.datasets import load_breast_cancer
|
||||
from sklearn.linear_model import LogisticRegression
|
||||
cancer = load_breast_cancer()
|
||||
import pandas as pd
|
||||
# Making a data frame
|
||||
cancerpd = pd.DataFrame(cancer.data, columns=cancer.feature_names)
|
||||
|
||||
fig, axes = plt.subplots(15,2,figsize=(10,20))
|
||||
malignant = cancer.data[cancer.target == 0]
|
||||
benign = cancer.data[cancer.target == 1]
|
||||
ax = axes.ravel()
|
||||
|
||||
for i in range(30):
|
||||
_, bins = np.histogram(cancer.data[:,i], bins =50)
|
||||
ax[i].hist(malignant[:,i], bins = bins, alpha = 0.5)
|
||||
ax[i].hist(benign[:,i], bins = bins, alpha = 0.5)
|
||||
ax[i].set_title(cancer.feature_names[i])
|
||||
ax[i].set_yticks(())
|
||||
ax[0].set_xlabel("Feature magnitude")
|
||||
ax[0].set_ylabel("Frequency")
|
||||
ax[0].legend(["Malignant", "Benign"], loc ="best")
|
||||
fig.tight_layout()
|
||||
plt.show()
|
||||
|
||||
import seaborn as sns
|
||||
correlation_matrix = cancerpd.corr().round(1)
|
||||
# use the heatmap function from seaborn to plot the correlation matrix
|
||||
# annot = True to print the values inside the square
|
||||
plt.figure(figsize=(15,8))
|
||||
sns.heatmap(data=correlation_matrix, annot=True)
|
||||
plt.show()
|
||||
|
||||
|
||||
!ec
|
||||
|
||||
!split
|
||||
===== Discussing the correlation data =====
|
||||
|
||||
In the above example we note two things. In the first plot we display
|
||||
the overlap of benign and malignant tumors as functions of the various
|
||||
features in the Wisconsing breast cancer data set. We see that for
|
||||
some of the features we can distinguish clearly the benign and
|
||||
malignant cases while for other features we cannot. This can point to
|
||||
us which features may be of greater interest when we wish to classify
|
||||
a benign or not benign tumour.
|
||||
|
||||
In the second figure we have computed the so-called correlation
|
||||
matrix, which in our case with thirty features becomes a $30\times 30$
|
||||
matrix.
|
||||
|
||||
We constructed this matrix using _pandas_ via the statements
|
||||
!bc pycod
|
||||
cancerpd = pd.DataFrame(cancer.data, columns=cancer.feature_names)
|
||||
!ec
|
||||
and then
|
||||
!bc pycod
|
||||
correlation_matrix = cancerpd.corr().round(1)
|
||||
!ec
|
||||
|
||||
Diagonalizing this matrix we can in turn say something about which
|
||||
features are of relevance and which are not. This leads us to
|
||||
the classical Principal Component Analysis (PCA) theorem with
|
||||
applications. This will be discussed later this semester ("week 43":"https://compphysics.github.io/MachineLearning/doc/pub/week43/html/week43-bs.html").
|
||||
|
||||
|
||||
|
||||
!split
|
||||
===== Other measures in classification studies: Cancer Data again =====
|
||||
!bc pycod
|
||||
import matplotlib.pyplot as plt
|
||||
import numpy as np
|
||||
from sklearn.model_selection import train_test_split
|
||||
from sklearn.datasets import load_breast_cancer
|
||||
from sklearn.linear_model import LogisticRegression
|
||||
|
||||
# Load the data
|
||||
cancer = load_breast_cancer()
|
||||
|
||||
X_train, X_test, y_train, y_test = train_test_split(cancer.data,cancer.target,random_state=0)
|
||||
print(X_train.shape)
|
||||
print(X_test.shape)
|
||||
# Logistic Regression
|
||||
logreg = LogisticRegression(solver='lbfgs')
|
||||
logreg.fit(X_train, y_train)
|
||||
|
||||
from sklearn.preprocessing import LabelEncoder
|
||||
from sklearn.model_selection import cross_validate
|
||||
#Cross validation
|
||||
accuracy = cross_validate(logreg,X_test,y_test,cv=10)['test_score']
|
||||
print(accuracy)
|
||||
print("Test set accuracy with Logistic Regression: {:.2f}".format(logreg.score(X_test,y_test)))
|
||||
|
||||
import scikitplot as skplt
|
||||
y_pred = logreg.predict(X_test)
|
||||
skplt.metrics.plot_confusion_matrix(y_test, y_pred, normalize=True)
|
||||
plt.show()
|
||||
y_probas = logreg.predict_proba(X_test)
|
||||
skplt.metrics.plot_roc(y_test, y_probas)
|
||||
plt.show()
|
||||
skplt.metrics.plot_cumulative_gain(y_test, y_probas)
|
||||
plt.show()
|
||||
|
||||
!ec
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
!split
|
||||
===== Gradient descent and Logistic regression =====
|
||||
|
||||
We complete these examples by adding a simple code for
|
||||
Logistic regression. Note the more general approach with a class for
|
||||
the method. Here we use a so-called _AND_ gate for our data set.
|
||||
|
||||
!bc pycod
|
||||
import numpy as np
|
||||
class LogisticRegression:
|
||||
def __init__(self, learning_rate=0.01, num_iterations=1000):
|
||||
self.learning_rate = learning_rate
|
||||
self.num_iterations = num_iterations
|
||||
self.beta_logreg = None
|
||||
def sigmoid(self, z):
|
||||
return 1 / (1 + np.exp(-z))
|
||||
def GDfit(self, X, y):
|
||||
n_data, num_features = X.shape
|
||||
self.beta_logreg = np.zeros(num_features)
|
||||
for _ in range(self.num_iterations):
|
||||
linear_model = X @ self.beta_logreg
|
||||
y_predicted = self.sigmoid(linear_model)
|
||||
# Gradient calculation
|
||||
gradient = (X.T @ (y_predicted - y))/n_data
|
||||
# Update beta_logreg
|
||||
self.beta_logreg -= self.learning_rate*gradient
|
||||
def predict(self, X):
|
||||
linear_model = X @ self.beta_logreg
|
||||
y_predicted = self.sigmoid(linear_model)
|
||||
return [1 if i >= 0.5 else 0 for i in y_predicted]
|
||||
# Example usage
|
||||
if __name__ == "__main__":
|
||||
# Sample data
|
||||
X = np.array([[0, 0], [1, 0], [0, 1], [1, 1]])
|
||||
y = np.array([0, 0, 0, 1]) # This is an AND gate
|
||||
model = LogisticRegression(learning_rate=0.01, num_iterations=1000)
|
||||
model.GDfit(X, y)
|
||||
predictions = model.predict(X)
|
||||
print("Predictions:", predictions)
|
||||
!ec
|
||||
|
||||
Executable
+3
@@ -0,0 +1,3 @@
|
||||
#!/bin/sh
|
||||
doconce clean
|
||||
rm -rf *.pdf *.tex ipynb*.tar.gz *.html ._*.html *~ reveal.js Trash README.txt
|
||||
Executable
+79
@@ -0,0 +1,79 @@
|
||||
#!/bin/sh
|
||||
set -x
|
||||
|
||||
function system {
|
||||
"$@"
|
||||
if [ $? -ne 0 ]; then
|
||||
echo "make.sh: unsuccessful command $@"
|
||||
echo "abort!"
|
||||
exit 1
|
||||
fi
|
||||
}
|
||||
|
||||
if [ $# -eq 0 ]; then
|
||||
echo 'bash make.sh slides1|slides2'
|
||||
exit 1
|
||||
fi
|
||||
|
||||
name=$1
|
||||
rm -f *.tar.gz
|
||||
|
||||
opt="--encoding=utf-8"
|
||||
# Note: Makefile examples contain constructions like ${PROG} which
|
||||
# looks like Mako constructions, but they are not. Use --no_mako
|
||||
# to turn off Mako processing.
|
||||
opt="--no_mako"
|
||||
|
||||
rm -f *.aux
|
||||
|
||||
|
||||
html=${name}-reveal
|
||||
system doconce format html $name --pygments_html_style=perldoc --keep_pygments_html_bg --html_links_in_new_window --html_output=$html $opt
|
||||
system doconce slides_html $html reveal --html_slide_theme=beige
|
||||
|
||||
# Plain HTML documents
|
||||
|
||||
html=${name}-solarized
|
||||
system doconce format html $name --pygments_html_style=perldoc --html_style=solarized3 --html_links_in_new_window --html_output=$html $opt
|
||||
system doconce split_html $html.html --method=space10
|
||||
|
||||
html=${name}
|
||||
system doconce format html $name --pygments_html_style=default --html_style=bloodish --html_links_in_new_window --html_output=$html $opt
|
||||
system doconce split_html $html.html --method=space10
|
||||
|
||||
# Bootstrap style
|
||||
html=${name}-bs
|
||||
system doconce format html $name --html_style=bootstrap --pygments_html_style=default --html_admon=bootstrap_panel --html_output=$html $opt
|
||||
system doconce split_html $html.html --method=split --pagination --nav_button=bottom
|
||||
|
||||
# IPython notebook
|
||||
system doconce format ipynb $name $opt
|
||||
|
||||
|
||||
# Publish
|
||||
dest=../../pub
|
||||
if [ ! -d $dest/$name ]; then
|
||||
mkdir $dest/$name
|
||||
mkdir $dest/$name/html
|
||||
mkdir $dest/$name/ipynb
|
||||
fi
|
||||
cp -r ${name}*.html ._${name}*.html reveal.js $dest/$name/html
|
||||
|
||||
# Figures: cannot just copy link, need to physically copy the files
|
||||
if [ -d fig-${name} ]; then
|
||||
if [ ! -d $dest/$name/html/fig-$name ]; then
|
||||
mkdir $dest/$name/html/fig-$name
|
||||
fi
|
||||
cp -r fig-${name}/* $dest/$name/html/fig-$name
|
||||
fi
|
||||
|
||||
cp ${name}.ipynb $dest/$name/ipynb
|
||||
ipynb_tarfile=ipynb-${name}-src.tar.gz
|
||||
if [ ! -f ${ipynb_tarfile} ]; then
|
||||
cat > README.txt <<EOF
|
||||
This IPython notebook ${name}.ipynb does not require any additional
|
||||
programs.
|
||||
EOF
|
||||
tar czf ${ipynb_tarfile} README.txt
|
||||
fi
|
||||
cp ${ipynb_tarfile} $dest/$name/ipynb
|
||||
Reference in New Issue
Block a user