OptimizationTechniques/Work 3/report/Work3_report.tex

%
% Optimization Techniques Work 3 report
%
% authors:
%   Χρήστος Χουτουρίδης ΑΕΜ 8997
%   cchoutou@ece.auth.gr


\documentclass[a4paper, 11pt]{AUTHReport}

% Document configuration
\AuthorName{Χρήστος Χουτουρίδης}
\AuthorAEM{8997}
\AuthorMail{cchoutou@ece.auth.gr}

%\CoAuthorName{CoAuthor Name}
%\CoAuthorAEM{AEM}
%\CoAuthorMail{CoAuthor Mail}

% \WorkGroup{Ομάδα Χ}

\DocTitle{3η Εργαστηριακή Άσκηση}
\DocSubTitle{Μέθοδος Μέγιστης Καθόδου με Προβολή}

\Department{Τμήμα ΗΜΜΥ. Τομέας Ηλεκτρονικής}
\ClassName{Τεχνικές Βελτιστοποίησης}

\InstructorName{Γ. Ροβιθάκης}
\InstructorMail{rovithak@auth.gr}

\CoInstructorName{Θ. Αφορόζη}
\CoInstructorMail{taforozi@ece.auth.gr}

\CurrentDate{\today}


\usepackage{enumitem}
\usepackage{tabularx}
\usepackage{array}
\usepackage{amssymb}
\usepackage{amsfonts}
\usepackage{amsmath}
\usepackage{commath}

\usepackage{float}
\usepackage[labelformat=empty]{subcaption}

\begin{document}

\InsertTitle

%\tableofcontents

\section{Εισαγωγή}
Η παρούσα εργασία αφορά το πρόβλημα της ελαχιστοποίησης μιας δοσμένης συνάρτησης πολλών μεταβλητών $f: \mathbb{R}^n \rightarrow \mathbb{R}$ με περιορισμούς χρησιμοποιώντας τη μέθοδο μέγιστης καθόδου με προβολή.
Η μέθοδος αυτή θα εκτελεστεί σε αντιπαραβολή με την αντίστοιχη μέθοδο χωρίς περιορισμούς από την προηγούμενη εργασία.
Για το λόγο αυτό χρησιμοποιούμε των κώδικα της προηγούμενης εργασίας με κάποιες τροποποιήσεις, όπως θα δούμε και παρακάτω.

\subsection{Παραδοτέα}
Τα παραδοτέα της εργασίας αποτελούνται από:
\begin{itemize}
    \item Την παρούσα αναφορά.
    \item Τον κατάλογο \textbf{scripts/}, που περιέχει τον κώδικα της MATLAB.
    \item Το \href{https://git.hoo2.net/hoo2/OptimizationTechniques/src/branch/master/Work%203}{σύνδεσμο} με το αποθετήριο που περιέχει όλο το project με τον κώδικα της MATLAB, της αναφοράς και τα παραδοτέα.
\end{itemize}

\subsection{Προγραμματιστική προσέγγιση}
Για τον προγραμματισμό και εκτέλεση των μεθόδων της παρούσας εργασίας έγινε χρήση της MATLAB.
Στον κατάλογο \textbf{scripts}, περιέχονται όλες οι μέθοδοι και οι τεχνικές υπολογισμού βημάτων με τη μορφή συναρτήσεων καθώς και scripts που τις καλούν.
Για κάθε ένα θέμα της εργασίας, υπάρχει το αντίστοιχο script που περιέχει τους υπολογισμούς, τις κλήσεις των μεθόδων και τη δημιουργία των διαγραμμάτων.
Για το πρώτο θέμα το αρχείο Script\_1\_SteepDesc.m για το δεύτερο το Script\_2\_SteepDesc\_Proj.m και ούτω καθεξής.
Η μέθοδος μέγιστης καθόδου (αρχείο: \textbf{method\_SteepDesc.m}) είναι η ίδια με αυτή της προηγούμενης εργασίας με τη μόνη διαφορά ότι τροποποιήθηκε ώστε η αντικειμενική συνάρτηση να δέχεται διάνυσμα ώς όρισμα και όχι δύο διαφορετικές μεταβλητές $x, y$.
Αυτό ακολουθήθηκε και για την έκδοση με προβολή και φυσικά είχε αντίκτυπο και στις υπόλοιπες συναρτήσεις, όπως η κλίση ή ο Εσσιανός.
Στην παρούσα εργασία η υλοποίηση του κώδικα ακολουθεί την προσέγγιση των προηγούμενων εργασιών.
Πιο συγκεκριμένα.

\subsection{Μέθοδοι επιλογής βήματος}
Εφόσον στην προηγούμενη εργασία οι μέθοδοι επιλογής βήματος ήταν ανεξάρτητες από την μέθοδο υπολογισμού του ελάχιστου και εφόσον χρησιμοποιούμε τον ίδιο κώδικά και στην παρούσα εργασία, αυτός ο τρόπος σχεδίασης παρέμεινε.
Ουσιαστικά για κάθε ένα τρόπο υπολογισμού του $\gamma_k$, υπάρχει αντίστοιχη συνάρτηση, με κοινό interface.
Αυτό έχει τη μορφή: \textit{\textbf{gamma\_<method>(f, grad\_f, dk, xk)}}, όπου το \textbf{f} είναι η αντικειμενική συνάρτηση, \textbf{grad\_f} η συνάρτηση κλίσης της, \textbf{dk} η τιμή της συνάρτησης κλίσης στο xk και \textbf{xk} το σημείο ενδιαφέροντος.
Έτσι οι μέθοδοι αντιγράφηκαν και εδώ για ολότητα, ακόμα και αν για την παρούσα εργασία χρησιμοποιείται μόνο το σταθερό βήμα $\gamma_k$.

%\subsection{Symbolic expression functions}
%Μία ακόμη προγραμματιστική τεχνική που ακολουθήθηκε είναι η χρήση \textbf{symbolic expression} για την αναπαράσταση της αντικειμενικής συνάρτησης.
%Η εξήγηση υπάρχει και στις προηγούμενες εργασίες αλλά την παραθέτουμε εδώ για ολότητα.
%Ο λόγος που επιλέχθηκε είναι η \textbf{δυνατότητα εξαγωγής ενός symbolic expression που αναπαριστά την κλίση $\nabla f$ και τον Εσσιανό $\nabla^2f$  μιας συνάρτησης} από την MATLAB, κάνοντας χρήση των εντολών \textit{gradient()} και \textit{hessian()}.
%Αν αντίθετα χρησιμοποιούσαμε απλές συναρτήσεις, πολυώνυμα ή lambdas για την αναπαράσταση των αντικειμενικών συναρτήσεων, τότε για τον υπολογισμό της κλίσης και του Εσσιανού θα έπρεπε:
%\begin{itemize}
%    \item Είτε να υπολογίζαμε αριθμητικά τις παραγώγους gradient και hessian μέσα στις μεθόδους, κάτι που θα εισήγαγε \textit{\textbf{αχρείαστο αριθμητικό σφάλμα}}.
%    \item Είτε να κάναμε χρήση δύο επιπλέων συναρτήσεων (ή πολυωνύμων) για την αναπαράσταση τους, κάτι που ουσιαστικά θα δημιουργούσε \textit{\textbf{πλεονασμό πληροφορίας εισόδου}} και άρα μεγαλύτερη πιθανότητα να κάνουμε λάθος.
%\end{itemize}
%Η αναπαράσταση όμως με χρήση symbolic expression είναι πιο “βαριά” όταν χρειάζεται να υπολογίσουμε την τιμή μιας συνάρτησης σε κάποιο σημείο (subs(expr, number)).
%Αυτό είναι κάτι που χρειάζεται εκτενώς στον κώδικά μας.
%Για το λόγο αυτό, ενώ η συνάρτηση δίνεται ως symbolic expression, μέσω αυτής υπολογίζονται αυτόματα η κλίση, ο Εσσιανός αλλά και οι “κανονικές” συναρτήσεις MATLAB που τις υλοποιούν.
%Έτσι έχουμε την ακριβή αναπαράσταση της κλίσης και του Εσσιανού ως συναρτήσεις χωρίς να πληρώνουμε το κόστος της subs().

\section{Απεικόνιση της συνάρτησης}
Η συνάρτηση με την οποία ασχολούμαστε στην παρούσα εργασία είναι η:
\boldmath
\begin{equation}
    f: \mathbb{R}^2 \rightarrow \mathbb{R}, f(x) = \frac{1}{3}{x_1}^2 + 3{x_2}^2, x = \begin{bmatrix} x_1 \\ x_2 \end{bmatrix}
\end{equation}
\label{eq:ObjectiveFunction}
Όπου $ -10 \leq x_1 \leq 5 $ και $ -8 \leq x_2 \leq 12 $ \\
\unboldmath
Στο παρακάτω σχήμα \ref{fig:plot3dFunction} φαίνεται η τρισδιάστατη απεικόνιση της συνάρτησης.
\InsertFigure{!h}{0.8}{fig:plot3dFunction}{../scripts/figures/Plot_Function.png}{Γραφική παράσταση της f}

Από το σχήμα μπορούμε πολύ εύκολα να διακρίνουμε ότι η συνάρτηση είναι κυρτή στο σύνολο των περιορισμών της εκφώνησης $ -10 \leq x_1 \leq 5 $ και $ -8 \leq x_2 \leq 12 $.
Για να πάρουμε μια καλύτερη αίσθηση για το που βρίσκεται το τοπικό ελάχιστο της $f$, παρακάτω παραθέτουμε ένα γράφημα με τις ισοβαρείς καμπύλες της $f$.
\InsertFigure{H}{0.8}{fig:plotContour}{../scripts/figures/Plot_Contour.png}{Ισοβαρείς της f}

Από το παραπάνω σχήμα \ref{fig:plotContour} φαίνονται και γραφικά οι μικρές κλίσης που παρουσιάζει η συνάρτηση κοντά στο ελάχιστο σημείο $(0,0)$.
Τα διαγράμματα για τη μέθοδο δημιουργούνται εκτελώντας το αρχείο \textbf{Script\_0\_Plots.m}

\section{Μέθοδος Μέγιστης Καθόδου χωρίς περιορισμούς - Θέμα 1}
Εφαρμόζοντας την μέθοδο μέγιστης καθόδου από την προηγούμενη εργασία, με ακρίβεια $\epsilon = 0.001$, για τα βήματα $\gamma_k$ της εκφώνησης, παρατηρούμε ότι η μέθοδος συγκλίνει στο ελάχιστο για μικρά $\gamma_k$ ενώ αποκλίνει για μεγάλα \boldmath$\gamma_k \geq 0.34$\unboldmath.
Από τις δοκιμές φαίνεται ότι το σημείο εκκίνησης δεν παίζει ρόλο και για αυτό επιλέξαμε να παραθέσουμε τα ευρήματά μας από το σημείο $(5,-5)$, για αντιπαραβολή με το επόμενο βήμα της εκφώνησης.
\InsertFigure{H}{0.6}{fig:StDes_Iter_o_gamma_2}{../scripts/figures/StDes_Iter_o_gamma_1.png}{Αριθμός επαναλήψεων για διαφορετικές τιμές $\gamma_k$ [Μέγιστη Κάθοδος].}
Επίσης παρατηρούμε ότι για μικρό \boldmath$\gamma_k = 0.1$ η σύγκλιση είναι ομαλή, ενώ για μεγάλο $\gamma_k = 0.3$ \unboldmath παρουσιάζει ταλάντωση κατά την σύγκλιση.
Παρακάτω στο σχήμα \ref{fig:StDes_gamma} παραθέτουμε την πορεία σύγκλισης και απόκλισης για τις διαφορετικές τιμές του $\gamma_k$.

\begin{figure}[ht]
    \centering
    % First row
    \begin{subfigure}{0.45\textwidth}
        \centering
        \includegraphics[width=\linewidth]{../scripts/figures/StDes_gamma_0.1.png}
        \caption{$\gamma_k = 0.1$}
        \label{fig:StDes_gamma_0.1}
    \end{subfigure}
    \hfill
    \begin{subfigure}{0.45\textwidth}
        \centering
        \includegraphics[width=\linewidth]{../scripts/figures/StDes_gamma_0.3.png}
        \caption{$\gamma_k = 0.3$}
        \label{fig:StDes_gamma_0.3}
    \end{subfigure}

    % Second row
    \vspace{1em}
    \begin{subfigure}{0.45\textwidth}
        \centering
        \includegraphics[width=\linewidth]{../scripts/figures/StDes_gamma_3.png}
        \caption{$\gamma_k = 3$}
        \label{fig:StDes_gamma_3}
    \end{subfigure}
    \hfill
    \begin{subfigure}{0.45\textwidth}
        \centering
        \includegraphics[width=\linewidth]{../scripts/figures/StDes_gamma_5.png}
        \caption{$\gamma_k = 5$}
        \label{fig:StDes_gamma_5}
    \end{subfigure}

    \caption{Σύγκριση της μεθόδου Steepest descent για διαφορετικά $\gamma_k$}
    \label{fig:StDes_gamma}
\end{figure}

\subsection{Μαθηματική ανάλυση}
Τα παραπάνω αποτελέσματα επιβεβαιώνονται και θεωρητικά.
Πιο συγκεκριμένα για τη σύγκλιση της μεθόδου μέγιστης καθόδου όπου το κάθε σημείο υπολογίζεται από την σχέση \boldmath$x_{k+1} = x_k - \gamma_k \nabla f(x_k)$, \unboldmath πρέπει να ισχύουν:
\begin{enumerate}
    \item H $f$ να είναι κυρτή.
    \item Η $f$ να είναι συνεχής και διαφορίσιμη και η κλίση της υπολογίσιμη.
    \item Για το βήμα υπολογισμού να ισχύει η σχέση:
    \begin{equation} 0 < \gamma_k < \frac{2}{L} \end{equation} \label{eq:gammaLimmit}
    Όπου $L$ το άνω φράγμα της Lipschitz για την κλίση $\nabla f(x)$ (αν είναι γνωστή), η οποία είναι η μέγιστη ιδιοτιμή του Εσσιανού και δίνεται από τη σχέση:
    \begin{equation}
        L = \max_{x} \{\lambda_{max} (H(x))\},  \quad H(x) = \begin{bmatrix} \dfrac{\theta^2 f(x)}{\theta x_i \theta x_j} \end{bmatrix}
    \end{equation}
    \label{eq:Lipschitz}
\end{enumerate}
\par
Έτσι για τη δική μας περίπτωση έχουμε:
Η κλίση της $f(x)$ είναι:
\[
\nabla f(x) = \begin{bmatrix} \frac{\partial f}{\partial x1} \\ \frac{\partial f}{\partial x2} \end{bmatrix} = \begin{bmatrix} \frac{2}{3}x_1 \\ 6x_2 \end{bmatrix}
\]
Ο Εσσιανός πίνακας $H(x)$ της $f(x)$ είναι:
\[
H(x) =
\begin{bmatrix}
    \frac{\partial^2 f}{\partial x_1^2} & \frac{\partial^2 f}{\partial x_1 \partial x_2} \\
    \frac{\partial^2 f}{\partial x_2 \partial x_1} & \frac{\partial^2 f}{\partial x_2^2}
\end{bmatrix}
=
\begin{bmatrix}
    \frac{2}{3} & 0 \\ 0 & 6
\end{bmatrix}
\]
Και εφόσον είναι διαγώνιος οι ιδιοτιμές του είναι οι τιμές της διαγωνίου. Δηλαδή:
\[
\det(H - \lambda I) = 0 \Leftrightarrow \\
\det(\begin{bmatrix}\frac{2}{3} - \lambda & 0 \\ 0 & 6 - \lambda \end{bmatrix}) = 0 \Leftrightarrow \\
\left(\frac{2}{3} - \lambda\right)(6 - \lambda) = 0 \Leftrightarrow \\
\lambda_{\min} = \frac{2}{3}, \quad \lambda_{\max} = 6.
\]
Έτσι από τις εξισώσεις \ref{eq:gammaLimmit} και \ref{eq:Lipschitz} προκύπτει τελικά:
\boldmath\[ 0 < \gamma_k < \frac{1}{3} \]\unboldmath
Βλέπουμε ότι από την ανάλυσή μας επιβεβαιώνουμε τις τιμές που βρήκαμε εμπειρικά από την εκτέλεση του αλγορίθμου.

\subsubsection{Εναλλακτικά}
Αν θέλαμε να βρούμε το κριτήριο σύγκλισης για το βήμα $\gamma_k$ ξεχωριστά για την κάθε διάσταση θα μπορούσαμε να θεωρήσουμε για 3ο κριτήριο πως για να συγκλίνει η μέθοδος θα πρέπει να ισχύει:
\[
    \norm{\frac{x_{k+1}}{x_k}} < 1
\]
Από την παραπάνω εξίσωση προκείπτει:


\section{Μέθοδος Μέγιστης Καθόδου με προβολή}
Πριν περάσουμε στις υπόλοιπες απαιτήσεις της εργασίας θα θέλαμε να παραθέσουμε κάποιες πληροφορίες για την υλοποίηση της μεθόδου μέγιστης καθόδου με προβολή (αρχείο: \textbf{method\_SteepDesc\_Proj.m}).
Η συνάρτηση αυτή δέχεται ως είσοδο την αντικειμενική συνάρτηση και την συνάρτηση κλίσης καθώς και το σημείο εκκίνησης $x_k$ και το βήμα $s_k$.
Με τη βοήθεια της συνάρτησης \textit{ProjectionPoint()} \textbf{παίρνει πρώτα την προβολή} του $x_k$ στο διάστημα των περιορισμών αν αυτό χρειάζεται \textbf{και έπειτα εφαρμόζει τον αλγόριθμο}.
Αυτό σημαίνει ότι μπορεί να χρησιμοποιηθεί και για σημεία εκκίνησης εκτός του συνόλου των περιορισμών.
Ο αλγόριθμος είναι παρόμοιος με αυτόν της προηγούμενης εργασίας με τη διαφορά ότι η διεύθυνση $d_k$ επιλέγεται από τη σχέση:
\[
    d_k = Pr_X\{ x_k - s_k \nabla f(x_k)\} - x_k
\]
Δηλαδή εφαρμόζουμε πρώτα τη μέθοδο μέγιστης καθόδου με βήμα $s_k$ στην κατεύθυνση $-\nabla f$ και έπειτα προβάλουμε το σημείο στο σύνολο $X$ και χρησιμοποιούμε αυτό ως διεύθυνση με βήμα $\gamma_k$.

\section{Μέθοδος Μέγιστης Καθόδου με προβολή $s_k = 5, \gamma_k = 0.5$ - Θέμα 2}
Εφαρμόζοντας τη μέθοδο για ακρίβεια $\epsilon = 0.01$, $s_k = 5$ και $\gamma_k = 0.5$ έχουμε:
\InsertFigure{H}{0.8}{fig:StDesProj_sk_5_gamma_0.5}{../scripts/figures/StDesProj_sk_5_gamma_0.5.png}{Μέθοδος μέγιστης καθόδου με προβολή για $s_k = 5, \gamma_k = 0.5$.}

Παρατηρούμε πως ενώ η μέθοδος ταλαντώνει και δεν συγκλίνει στο ελάχιστο \textbf{όπως και η αντίστοιχη εκτέλεση της μέγιστης καθόδου χωρίς περιορισμούς με το ίδιο} \boldmath$\gamma_k$\unboldmath.
Παρόλα αυτά όμως, η ταλάντωση λαμβάνει χώρα \textbf{εντός του συνόλου των περιορισμών} της εκφώνησης.

\section{Μέθοδος Μέγιστης Καθόδου με προβολή $s_k = 15, \gamma_k = 0.1$ - Θέμα 3}
Εφαρμόζοντας τη μέθοδο για ακρίβεια $\epsilon = 0.01$, $s_k = 15$ και $\gamma_k = 0.1$ έχουμε:
\InsertFigure{H}{0.8}{fig:StDesProj_sk_15_gamma_0.1}{../scripts/figures/StDesProj_sk_15_gamma_0.1.png}{Μέθοδος μέγιστης καθόδου με προβολή για $s_k = 15, \gamma_k = 0.1$.}

Εδώ παρατηρούμε πως ενώ το $\gamma_k$ έχει επιλεγεί στο εύρος που οδηγεί σε σύγκλιση, το αντίστοιχο βήμα $s_k$ είναι πολύ μεγάλο, με αποτέλεσμα η μέθοδος να ταλαντώνει και πάλι.
Αυτή τη φορά μόνο στον άξονα $x_2$.

\section{Μέθοδος Μέγιστης Καθόδου με προβολή $s_k = 0.1, \gamma_k = 0.2$ - Θέμα 4}
Αρχικά παρατηρούμε πως το σημείο δεν είναι εφικτό, καθώς είναι εκτός του συνόλου των περιορισμών της εκφώνησης.
Αυτό βέβαια δεν μας αποτρέπει από την εφαρμογή της μεθόδου, καθώς αρχικά μπορούμε να προβάλουμε το σημείο στο σύνολο και να εφαρμόσουμε τη μέθοδο έπειτα.
Ακόμα, αυτή τη φορά οι τιμές των βημάτων $s_k, \gamma_k$, έχουν επιλεγεί μέσα στο εύρος για το οποίο έχουμε σύγκλιση, επομένως αναμένουμε η μέθοδος να συγκλίνει στο ελάχιστο.
Εφαρμόζοντας τη μέθοδο για ακρίβεια $\epsilon = 0.01$ έχουμε:
\InsertFigure{H}{0.8}{fig:StDesProj_sk_0.1_gamma_0.2}{../scripts/figures/StDesProj_sk_0.1_gamma_0.2.png}{Μέθοδος μέγιστης καθόδου με προβολή για $s_k = 0.1, \gamma_k = 0.2$.}

\section{Συμπεράσματα}
...

\end{document}