WIP: Pre-final version

10_Introduction.tex

@@ -1,5 +1,7 @@
 \section{Introduction}
 
+\subsection{Motivation}
+
 Buildings are a major consumer of energy, with more than 25\% of the total
 energy consumed in the EU coming from residential
 buildings~\cite{tsemekiditzeiranakiAnalysisEUResidential2019}. Combined with a
@@ -18,10 +20,88 @@ Gaussian Processes have been previously used to model building dynamics, but
 they are usually limited by a fixed computational budget. This limits the
 approaches that can be taken for identification and update of said models.
 Learning \acrshort{gp} models have also been previously used in the context of
-autonomous racing cars, but there the Sparse \acrshort{gp} model was built on
-top of a white-box model and only responsible for fitting the unmodeled
-dynamics. 
+autonomous racing cars \cite{kabzanLearningBasedModelPredictive2019}, but there
+the Sparse \acrshort{gp} model was built on top of a white-box model and only
+responsible for fitting the unmodeled dynamics. 
 
-This project means to provide a further expansion of the use of black-box
-\acrlong{gp} Models in the context of building control, through online learning
-of building dynamics at new operating points as more data gets collected.
+\subsection{Previous Research}
+With the increase in computational power and availability of data  over time,
+the accessibility of data-driven methods for system identification and control
+has also risen significantly. 
+
+The idea of using Gaussian Processes as regression models for control of dynamic
+systems is not new, and has already been explored a number of times. A general
+description of their use, along with the necessary theory and some example
+implementations is given in~\cite{kocijanModellingControlDynamic2016}.
+In~\cite{pleweSupervisoryModelPredictive2020}, a \acrlong{gp} Model with a
+\acrlong{rq} Kernel is used for temperature set point optimization.
+
+Gaussian Processes for building control have also been studied before in the
+context of Demand Response~\cite{nghiemDatadrivenDemandResponse2017,
+jainLearningControlUsing2018}, where the buildings are used for their heat
+capacity in order to reduce the stress on energy providers during peak load
+times.
+
+There are, however, multiple limitations with these approaches. 
+In~\cite{nghiemDatadrivenDemandResponse2017} the model is only identified once,
+ignoring changes in weather or plant parameters, which could lead to different
+dynamics. This is addressed in \cite{jainLearningControlUsing2018} by
+re-identifying the model every two weeks using new information. Another
+limitation is that of the scalability of the \acrshort{gp}s, which become
+prohibitively expensive from a computational point of view when too much data is
+added.
+
+Outside of the context of building control, Sparse \acrlong{gp}es have been used
+in autonomous racing in order to complement the physics-based model by fitting
+the unmodeled dynamics of the
+system~\cite{kabzanLearningBasedModelPredictive2019}.
+
+\subsection{Contribution}
+
+The ability to learn the plant's behaviour in new regions is very helpful in
+maintaining model performance over time, as its behaviour starts deviating and
+the original identified model goes further and further into the extrapolated
+regions.
+
+This project tries to combine the use of online learning control schemes with
+\acrlong{gp} Models through implementing \acrlong{svgp} Models. \acrshort{svgp}s
+provide means of extending the use of \acrshort{gp}s to larger datasets, thus
+enabling the periodic re-training of the model to include all the historically
+available data.
+
+\subsection{Project Outline}
+
+The main body of work can be divided in two parts: the development of a computer
+model of the \pdome\ building and the synthesis, validation and comparison of
+multiple control schemes using both classical \acrshort{gp}s, as well as
+\acrshort{svgp}s.
+
+Section~\ref{sec:gaussian_processes} provides the mathematical background for
+understanding \acrshort{gp}s, as well as the definition in very broad strokes of
+\acrshort{svgp}s and their differences from the classical implementation of
+\acrlong{gp}es. This information is later used for comparing their performances
+and outlining their respective pros and cons.
+
+Section~\ref{sec:CARNOT} goes into the details of the implementation of the
+\pdome\ computer model. The structure of the CARNOT model is described, and all
+the parameters required for the simulation are either directly sourced from
+existing literature, computed from available values or estimated using publicly
+available data from \textit{Google Maps}~\cite{GoogleMaps}.
+
+Following that, Section~\ref{sec:hyperparameters} discusses the choice of the
+hyperparameters for the \acrshort{gp} and \acrshort{svgp} models respectively in
+the context of their multi-step ahead simulation performance.
+
+The \acrlong{ocp} implemented in the \acrshort{mpc} controller is described in
+Section~\ref{sec:mpc_problem}, followed by a description of the complete
+controller implementation in Python in Section~\ref{sec:implementation}.
+
+Section~\ref{sec:results} presents and analyzes the results of full-year
+simulations for both \acrshort{gp} and \acrshort{svgp} models. A few variations
+on the initial \acrshort{svgp} model are further analyzed in order to identify
+the most important parameteres for long time operation.
+
+Finally, Section~\ref{sec:further_research} discusses the possible improvements
+to both the current \acrshort{gp} and \acrshort{svgp} implementations.
+
+\clearpage

20_Previous_Research.tex

@@ -1,44 +0,0 @@
-\section{Previous Research}
-With the increase in computational power and availability of data  over time,
-the accessibility of data-driven methods for System Identification and Control
-has also risen significantly. 
-
-The idea of using Gaussian Processes as regression models for control of dynamic
-systems is not new, and has already been explored a number of times. A general
-description of their use, along with the necessary theory and some example
-implementations is given in~\cite{kocijanModellingControlDynamic2016}.
-In~\cite{pleweSupervisoryModelPredictive2020} a \acrlong{gp} Model with a
-\acrlong{rq} Kernel is used for temperature set point optimization.
-
-Gaussian Processes for building control have also been studied before in the
-context of Demand Response~\cite{nghiemDatadrivenDemandResponse2017,
-jainLearningControlUsing2018}, where the buildings are used for their heat
-capacity in order to reduce the stress on energy providers during peak load
-times.
-
-There are, however, multiple limitations with these approaches. 
-In~\cite{nghiemDatadrivenDemandResponse2017} the model is only identified once,
-ignoring changes in weather or plant parameters, which could lead to different
-dynamics. This is addressed in \cite{jainLearningControlUsing2018} by
-re-identifying the model every two weeks using new information. Another
-limitation is that of the scalability of the \acrshort{gp}s, which become
-prohibitively expensive from a computational point of view when too much data is
-added.
-
-Outside of the context of building control, Sparse \acrlong{gp}es have been used
-in autonomous racing in order to complement the physics-based model by fitting
-the unmodeled dynamics of the
-system~\cite{kabzanLearningBasedModelPredictive2019}.
-
-The ability to learn the plant's behaviour in new regions is very helpful in
-maintaining model performance over time, as its behaviour starts deviating and
-the original identified model goes further and further into the extrapolated
-regions.
-
-
-This project will therefore try to combine the use of online learning schemes
-with \acrlong{gp}es by implementing \acrlong{svgp}es, which provide means of
-employing \acrshort{gp} Models on larger datasets, and re-training the models
-every day at midnight to include all the historically available data.
-
-\clearpage

30_Gaussian_Processes_Background.tex

@@ -36,7 +36,7 @@ to the given prior:
     \begin{bmatrix}
         \mathbf{f} \\
         \mathbf{f_*} \\
-    \end{bmatrix} = 
+    \end{bmatrix} \sim
     \mathcal{N}\left(
         \mathbf{0}, 
         \begin{bmatrix}
@@ -53,7 +53,7 @@ In the case of noisy observations, assuming $y = f + \epsilon$  with $\epsilon
     \begin{bmatrix}
         \mathbf{y} \\
         \mathbf{f_*} \\
-    \end{bmatrix} = 
+    \end{bmatrix} \sim
     \mathcal{N}\left(
         \mathbf{0}, 
         \begin{bmatrix}
@@ -69,7 +69,7 @@ which, for the rest of the section, will be used in the abbreviated form:
     \begin{bmatrix}
         \mathbf{y} \\
         \mathbf{f_*} \\
-    \end{bmatrix} = 
+    \end{bmatrix} \sim
     \mathcal{N}\left(
         \mathbf{0}, 
         \begin{bmatrix}
@@ -96,7 +96,7 @@ value that gets maximized is the log of Equation~\ref{eq:gp_likelihood}, the log
 marginal likelihood:
 
 \begin{equation}\label{eq:gp_log_likelihood}
-    log(p(y)) = - \frac{1}{2}\log{\left(
+    \log(p(y)) = - \frac{1}{2}\log{\left(
                                 \det{\left(
                                         K + \sigma_n^2I
                                 \right)}
@@ -172,7 +172,7 @@ $\mathbf{x'}$'s dimensions:
 where $w_d = \frac{1}{l_d^2}; d = 1 ,\dots, D$, with D being the dimension of the
 data.
 
-The special case of $\Lambda^{-1} = \text{diag}{\left([l_1^{-2},\dots,l_D^{-2}]\right)}$
+This special case of $\Lambda^{-1} = \text{diag}{\left([l_1^{-2},\dots,l_D^{-2}]\right)}$
 is equivalent to implementing different lengthscales on different regressors.
 This can be used to asses the relative importance of each regressor through the
 value of the hyperparameters. This is the \acrfull{ard} property.
@@ -207,7 +207,7 @@ without inquiring the penalty of inverting the covariance matrix. An overview
 and comparison of multiple methods is given
 at~\cite{liuUnderstandingComparingScalable2019}.
 
-For the scope of this project the choice of using the \acrfull{svgp} models has
+For the scope of this project, the choice of using the \acrfull{svgp} models has
 been made, since it provides a very good balance of scalability, capability,
 robustness and controllability~\cite{liuUnderstandingComparingScalable2019}.
 
@@ -230,9 +230,9 @@ $f(X_s)$, usually denoted as $f_s$ is introduced, with the requirement that this
 new dataset has size $n_s$  smaller than the size $n$ of the original dataset. 
 
 The $X_s$ are called \textit{inducing locations}, and $f_s$ --- \textit{inducing
-random variables}. They are said to summarize the data in the sense that a model
-trained on this new dataset should be able to generate the original dataset with
-a high probability.
+random variables}. They summarize the data in the sense that a model trained on
+this new dataset should be able to generate the original dataset with a high
+probability.
 
 The multivariate Gaussian distribution is used to establish the relationship
 between $f_s$ and $f$, which will serve the role of the prior, now called the
@@ -242,7 +242,7 @@ sparse prior:
     \begin{bmatrix}
         \mathbf{f}(X) \\
         \mathbf{f}(X_s) \\
-    \end{bmatrix} = 
+    \end{bmatrix} \sim
     \mathcal{N}\left(
         \mathbf{0}, 
         \begin{bmatrix}
@@ -267,8 +267,8 @@ computationally tractable on larger sets of data.
 The following derivation of the \acrshort{elbo} is based on the one presented
 in~\cite{yangUnderstandingVariationalLower}.
 
-Assume $X$ to be the observations, and $Z$ the set of hidden (latent)
-variables --- the parameters of the \acrshort{gp} model. The posterior
+Assume $X$ to be the observations, and $Z$ the set parameters of the
+\acrshort{gp} model, also known as the latent variables. The posterior
 distribution of the hidden variables can be written as follows:
 
 \begin{equation}
@@ -337,7 +337,10 @@ The \acrfull{noe} uses the past predictions $\hat{y}$ for future predictions:
     f(w(k-1),\dots,w(k-l_w),\hat{y}(k-1),\dots,\hat{y}(k-l_y),u(k-1),\dots,u(k-l_u))
 \end{equation}
 
-The \acrshort{noe} structure is therefore a \textit{simulation model}.
+Due to its use for multi-step ahead simulation of system behaviour, as opposed
+to only predicting one state ahead using current information, the \acrshort{noe}
+structure can be considered a \textit{simulation model}.
+
 
 In order to get the best simulation results from a \acrshort{gp} model, the
 \acrshort{noe} structure would have to be employed. Due to the high algorithmic

40_CARNOT_model.tex

@@ -94,7 +94,7 @@ Table~\ref{tab:GIMP_measurements}:
         \hline
         Object & Size [px] & Size[mm] & Size[m]\\
         \hline \hline
-        acrshort{hvac} height & 70 & 2100 & 2.1 \\
+        \acrshort{hvac} height & 70 & 2100 & 2.1 \\
         Building height & 230 & 6900 & 6.9 \\
         Stem wall & 45 & 1350 & 1.35 \\
         Dome height & 185 & 5550 & 5.55 \\
@@ -174,13 +174,13 @@ The stem of the \pdome\ is approximated to a cube of edge 25m, and its volume ca
 therefore be calculated as: 
 
 \begin{equation}
-    V_s = l_s^2 * h_s
+    V_s = l_s^2 \times h_s
 \end{equation}
 
 The total volume of the building is then given as: 
 
 \begin{equation}
-    V = V_d + V_s = \frac{1}{6} \pi h (3r^2 + h^2) + l_s^2 * h_s
+    V = V_d + V_s = \frac{1}{6} \pi h (3r^2 + h^2) + l_s^2 \times h_s
 \end{equation}
 
 Numerically, considering a dome diameter of 28m, a dome height of 5.55m and a stem
@@ -199,7 +199,7 @@ Exchange Rate, presented in Section~\ref{sec:Air_Exchange_Rate}.
 The main function of the \pdome\ building is serving as a classroom for around
 one hundred students. It has wood furniture consisting of chairs and tables, as
 well as a wooden stage in the center of the building, meant to be used for
-presentations. The building also contains a smaller room, housing the all the
+presentations. The building also contains a smaller room, housing all the
 necessary technical equipment (cf. Figure~\ref{fig:Google_Maps_Streetview}).
 
 The most accurate way of including information on the furniture in the CARNOT
@@ -233,7 +233,15 @@ has a surface area of:
 
 \begin{equation}
     S_f = \frac{1}{4} \cdot 1.8 \left[\frac{\text{m}^2}{\text{m}^2}\right]
-    \cdot 625\ \left[\text{m}^2\right] = 140\ \left[\text{m}^2\right]
+    \cdot 625\ \left[\text{m}^2\right] = 281.25 \left[\text{m}^2\right]
+\end{equation}
+
+Since a wall has two large faces, the total surface gets divided between them
+for a surface of each face of: 
+
+\begin{equation}
+    S_{f, \text{face}} = 
+    \frac{1}{2} \times S_f \approx 140 \left[\text{m}^2\right]
 \end{equation}
 
 \subsubsection*{Mass}
@@ -270,8 +278,9 @@ The last parameter of the furniture equivalent wall is computed by dividing its
 volume by the surface:
 
 \begin{equation}
-    h_f = \frac{V_f}{S_f} = \frac{10.41}{140} \left[\frac{m^3}{m^2}\right] =
-    0.075\ \left[\text{m}\right] = 7.5\ \left[\text{cm}\right]
+    h_f = \frac{V_f}{S_{f, \text{face}}} 
+    = \frac{10.41}{140} \left[\frac{m^3}{m^2}\right]
+    = 0.075\ \left[\text{m}\right] = 7.5\ \left[\text{cm}\right]
 \end{equation}
 
 
@@ -284,7 +293,7 @@ parameters for each of the CARNOT nodes' properties.
 
 \subsubsection{Windows}
 
-The windows are supposed to be made of two glass panes of thickness 4mm each.
+The windows are made of two glass panes of thickness 4mm each.
 
 The values of the heat transfer coefficient (U-factor) can vary greatly for
 different window technologies, but an educated guess can be made on the lower
@@ -345,8 +354,8 @@ Table~\ref{tab:material_properties}:
 \centering
     \begin{tabular}{||c c c c||}
         \hline
-        Material & Thermal Conductivity $[k]$ $[\frac{W}{mK}]$ & Specific Heat
-        Capacity $[c]$ $[\frac{J}{kgK}]$ & Density $[\rho]$ $[\frac{kg}{m^3}]$
+        Material & Thermal Conductivity $k$ $[\frac{W}{mK}]$ & Specific Heat
+        Capacity $c$ $[\frac{J}{kgK}]$ & Density $\rho$ $[\frac{kg}{m^3}]$
         \\
         \hline \hline
         Plywood & 0.12 & 1210 & 540 \\
@@ -375,7 +384,7 @@ the cases where the full cooling capacity is not required.
 \subsubsection*{Ventilation}
 
 According to the manufacturer manual \cite{aermecRoofTopManuelSelection}, the
-\acrshort{hvac} unit's external fan has an air debit ranging between 4900
+\acrshort{hvac} unit's external fan has an air flow ranging between 4900
 $\text{m}^3/\text{h}$ and 7000 $\text{m}^3/\text{h}$.
 
 \subsubsection*{Air Exchange Rate}\label{sec:Air_Exchange_Rate}
@@ -467,12 +476,13 @@ compressors. The instruction manual of the
 \acrshort{hvac}~\cite{aermecRoofTopManuelSelection} notes that in summer only
 one of the compressors is running. This allows for a larger \acrshort{eer} value
 and thus better performance. We can see that this is the case for most of the
-experiment, where the power consumption caps at around 6~kW. There are, however,
-moments during the first part of the experiment where the power momentarily
-peaks over the 6~kW limit, and goes as high as around 9~kW. This most probably
-happens when the \acrshort{hvac} decides that the difference between the set
-point temperature and the actual measured values is too large to compensate with
-a single compressor.
+experiment, where the power consumption caps at around 6~kW, which is the
+maximum power consumption of one compressor. There are, however, moments during
+the first part of the experiment where the power momentarily peaks over the 6~kW
+limit, and goes as high as around 9~kW. This most probably happens when the
+\acrshort{hvac} decides that the difference between the set point temperature
+and the actual measured values is too large to compensate with a single
+compressor.
 
 Figure~\ref{fig:Polydome_exp7_settemp} presents the values of the set point
 temperature and the measured internal temperature. 
@@ -506,7 +516,7 @@ pressure, wind speed and direction, etc.  A detailed overview of each
 measurement necessary for a simulation is given in the CARNOT user
 manual~\cite{CARNOTManual}.
 
-In order to compare the CARNOT model's performance to that of the real \pdome\,
+In order to compare the CARNOT model's performance to that of the real \pdome,
 it is necessary to simulate the CARNOT model under the same set of conditions as
 the ones present during the experimental data collection. In order to do this,
 all the missing values that are required by the simulation have to be filled. In
@@ -540,12 +550,12 @@ All the other parameters related to solar irradiance, such as the in-plane
 irradiance components, in-plane diffuse irradiance from the sky and the ground
 are computed using the Python pvlib.
 
-The values that cannot be either calculated or approximated from the available
+The values that can be neither calculated nor approximated from the available
 data, such as the precipitation, wind direction, incidence angles in place of
 vertical and main/secondary surface axis, have been replaced with the default
 CARNOT placeholder value of -9999. The relative humidity, cloud index, pressure
-and wind speed have been fixed to 50\%, 0.5, 96300 Pa, 0 $\text{m}/\text{s}$
-respectively.
+and wind speed have been fixed to the default CARNOT values of 50\%, 0.5, 96300
+Pa, 0 $\text{m}/\text{s}$ respectively.
 
 \subsubsection{\pdome\ and CARNOT model comparison}\label{sec:CARNOT_comparison}
 

50_Choice_of_Hyperparameters.tex

@@ -1,4 +1,4 @@
-\section{Choice of Hyperparameters}
+\section{Choice of Hyperparameters}~\label{sec:hyperparameters}
 
 This section will discuss and try to validate the choice of all the
 hyperparameters necessary for the training of a \acrshort{gp} model to capture
@@ -9,16 +9,14 @@ behaviour directly from data. This comes in contrast to white-box and grey-box
 modelling techniques, which require much more physical insight into the plant's
 behaviour.
 
-The advantage of black-box models lies in the lack of physical parameters to be
-fitted. On the flip side, this versatility of being able to fit much more
-complex models purely on data comes at the cost of having to properly define the
-model hyperparameters: the number of regressors, the number of autoregressive
-lags for each class of inputs, the shape of the covariance function have to be
-taken into account when designing a \acrshort{gp} model. These choices have
-direct influence on the resulting model behaviour and where it can be
-generalized, as well as indirect influence in the form of more time consuming
-computations in the case of larger number of regressors and more complex kernel
-functions.
+On the flip side, this versatility comes at the cost of having to properly
+define the model hyperparameters: the number of regressors, the number of
+autoregressive lags for each class of inputs, the shape of the covariance
+function have to be taken into account when designing a \acrshort{gp} model.
+These choices have direct influence on the resulting model behaviour and where
+it can be generalized, as well as indirect influence in the form of more time
+consuming computations in the case of larger number of regressors and more
+complex kernel functions.
 
 As described in Section~\ref{sec:gp_dynamical_system}, for the purpose of this
 project, the \acrlong{gp} model will be trained using the \acrshort{narx}
@@ -64,7 +62,12 @@ always benefit the \acrshort{gp} model at least not comparably to the
 additional computational complexity.
 
 For the exogenous inputs the choice has therefore been made to take the
-\textit{Global Solar Irradiance} and \textit{Outside Temperature Measurement}.
+\textit{Global Solar Irradiance}\footnotemark and \textit{Outside Temperature
+Measurement}.
+
+\footnotetext{Using the \acrshort{ghi} makes sense in the case of the \pdome\
+building which has windows facing all directions, as opposed to a room which
+would have windows on only one side.}
 
 \subsection{The Kernel}
 
@@ -104,7 +107,7 @@ three lengthscales apart.
     \begin{tabular}{||c c ||}
         \hline
         $\norm{\mathbf{x} - \mathbf{x}'}$ &
-        $\exp{(-\frac{1}{2}*\frac{\norm{\mathbf{x} - \mathbf{x}'}^2}{l^2})}$ \\
+        $\exp{(-\frac{1}{2}\times\frac{\norm{\mathbf{x} - \mathbf{x}'}^2}{l^2})}$ \\
         \hline \hline
         $1l$ & 0.606 \\
         $2l$ & 0.135 \\
@@ -133,6 +136,8 @@ and the variance for different combinations of the exogenous input lags ($l_w$),
 the controlled input lags ($l_u$) and the output lags ($l_y$) for a classical
 \acrshort{gp} model. 
 
+% TODO: [Lengthscales] Explain lags and w1,1, w2, etc.
+
 \begin{table}[ht]
 %\vspace{-8pt}
 \centering
@@ -165,15 +170,15 @@ the controlled input lags ($l_u$) and the output lags ($l_y$) for a classical
 \label{tab:GP_hyperparameters}
 \end{table}
 
-In general, the results of Table~\ref{tab:GP_hyperparameters} show that the
-past outputs are important when predicting future values. Of importance is also
-the past inputs, with the exception of the models with very high variance, where
-the relative importances stay almost constant across all the inputs. For the
-exogenous inputs, the outside temperature ($w2$) is generally more important
-than the solar irradiation ($w1$). In the case of more autoregressive lags for
-the exogenous inputs the more recent information is usually more important in
-the case of the solar irradiation, while the second-to-last measurement is
-preffered for the outside temperature.
+In general, the results of Table~\ref{tab:GP_hyperparameters} show that the past
+outputs are important when predicting future values. Of importance is also the
+past inputs, with the exception of the models with very high variance, where the
+relative importances stay almost constant across all the inputs. For example, it
+can be seen that for the exogenous inputs, the outside temperature ($w2$) is
+generally more important than the solar irradiation ($w1$). In the case of more
+autoregressive lags for the exogenous inputs the more recent information is
+usually more important in the case of the solar irradiation, while the
+second-to-last measurement is preffered for the outside temperature.
 
 For the classical \acrshort{gp} model the appropriate choice of lags would be
 $l_u = 1$ and $l_y = 3$ with $l_w$ taking the values of either 1, 2 or 3,
@@ -348,8 +353,6 @@ metrics, as well as being cheaper from a computational perspective. In order to
 make a more informed choice for the best hyperparameters, the simulation
 performance of all three combinations has been analysed.
 
-\clearpage
-
 \begin{table}[ht]
 %\vspace{-8pt}
 \centering
@@ -381,7 +384,6 @@ other combinations scoring much worse on the \acrshort{msll} and \acrshort{lpd}
 loss functions. This has, therefore, been chosen as the model for the full year
 simulations.
 
-
 \subsection{Validation of hyperparameters}\label{sec:validation_hyperparameters}
 
 The validation step has the purpose of testing the viability of the trained
@@ -449,10 +451,13 @@ this proves to be the best simulation model.
 
 Lastly, \model{3}{1}{3} has a much worse simulation performance than the other
 two models. This could hint at an over fitting of the model on the training data.
+
+\clearpage
+
 This is consistent with the results found in Table~\ref{tab:GP_loss_functions}
 for the \acrshort{rmse} and \acrshort{smse}, as well as can be seen in
 Appendix~\ref{apx:hyperparams_gp}, Figure~\ref{fig:GP_313_test_validation},
-where the model has much worse performance on the testing dataset predictions
+where \model{3}{1}{3} has much worse performance on the testing dataset predictions
 than the other two models.
 
 The performance of the three models in simulation mode is consistent with the
@@ -471,8 +476,6 @@ while still having good performance on all loss functions. In implementation,
 however, this model turned out to be very unstable, and the more conservative
 \model{1}{1}{3} model was used instead.
 
-\clearpage
-
 \subsubsection{Sparse and Variational Gaussian Process}
 
 For the \acrshort{svgp} models, only the performance of \model{1}{2}{3} was

60_The_MPC_Problem.tex

@@ -1,6 +1,6 @@
 \section{The MPC Problem}\label{sec:mpc_problem}
 
-The \acrlong{ocp} to be solved was chosen in such a way as to make
+The \acrlong{ocp} to be solved was chosen in such a way as to make the
 analysis of the models' performances more straightforward. The objective is
 tracking a defined reference temperature as close as possible, while ensuring
 the heat input stays within the HVAC capacity. The \textit{zero-variance} method
@@ -42,7 +42,8 @@ buildings based on the rolling 48h average outside temperature.
 \begin{figure}[ht]
     \centering
     \includegraphics[width = \textwidth]{Images/sia_180_2014.png}
-    \caption{The SIA 180:2014 norm for residential building temperatures}
+    \caption{The SIA 180:2014 norm for residential building
+    temperatures~\cite{sia180:2014ProtectionThermiqueProtection2014}}
     \label{fig:sia_temperature_norm}
 \end{figure}
 

70_Implementation.tex

@@ -21,13 +21,13 @@ in the \acrshort{wdb} object is given in Section~\ref{sec:CARNOT_WDB}.
 \subsection{Simulink Model}
 
 The secondary functions of the Simulink model is the weather prediction, as well
-as communication with the Python controller. A complete schema of the Simulink
-setup is presented in Figure~\ref{fig:Simulink_complete}.
+as communication with the Python controller. A complete schematic of the
+Simulink setup is presented in Figure~\ref{fig:Simulink_complete}.
 
 \begin{figure}[ht]
     \centering
     \includegraphics[width = 0.75\textwidth]{Images/polydome_python.pdf}
-    \caption{Simulink Schema of the Complete Simulation}
+    \caption{Simulink diagram of the Complete Simulation}
     \label{fig:Simulink_complete}
 \end{figure}
 
@@ -158,7 +158,7 @@ Let $w_l$, $u_l$, and $y_l$ be the lengths of the state vector components
 $\mathbf{w}$, $\mathbf{u}$, $\mathbf{y}$ (cf. Equation~\ref{eq:components}).
 Also, let X be the matrix of all the system states over the optimization horizon
 and W be the matrix of the predicted disturbances for all the future steps. The
-original \acrlong{ocp} can be rewritten as:
+original \acrlong{ocp} can be rewritten using index notation as:
 
 \begin{subequations}\label{eq:sparse_optimal_control_problem}
     \begin{align}
@@ -174,7 +174,7 @@ original \acrlong{ocp} can be rewritten as:
     \end{align}
 \end{subequations}
 
-\subsection{Python server}
+\subsection{Model Identification and Update using a Python server}
 
 The Python server is responsible for the control part of the simulation. It
 delegates which controller is active, is responsible for training and updating
@@ -182,16 +182,16 @@ the \acrshort{gp} and \acrshort{svgp} models, as well as keeping track of all
 the intermediate results for analysis.
 
 In the beginning of the experiment there is no information available on the
-building's thermal behaviour. For this part of the simulation, the controller
+building's thermal behaviour. For this part of the simulation, the server
 switches to a \acrshort{pi} controller with random disturbances until it gathers
-enough data to train a \acrshort{gp} model. This ensured that the building is
+enough data to train a \acrshort{gp} model. This ensures that the building is
 sufficiently excited to capture its dynamics, while maintaining the temperature
 within an acceptable range (~15 --- 25 $\degree$C).
 
 Once enough data has been captured all the features are first scaled to the
 range [-1, 1] in order to reduce the possibility of numerical instabilities. The
-Python controller then trains the \acrshort{gp} model and switches to tracking
-the appropriate SIA 180:2014 reference temperature (cf.
+Python server then trains the \acrshort{gp} model and switches to tracking the
+appropriate SIA 180:2014 reference temperature (cf.
 Section~\ref{sec:reference_temperature}).
 
 For the case of the \acrshort{svgp}, a new model is trained once enough data is

80_Results.tex

@@ -1,17 +1,17 @@
-\section{Results}\label{sec:results}
+\section{Full-year Simulation Results}\label{sec:results}
 
 
 This section focuses on the presentation and interpretation of the year-long
-simulation of the control schemes present previously. All the control schemes
-analysed in this Section have been done with a sampling time of 15 minutes and a
-control horizon of 8 steps.
+simulation of the control schemes presented previously. All the control schemes
+analysed in this Section have used a sampling time of 15 minutes and a control
+horizon of 8 steps.
 
 Section~\ref{sec:GP_results} analyses the results of a conventional
-\acrlong{gp} Model trained on the first five days of gathered data. The models
+\acrlong{gp} Model trained on the first five days of gathered data. The model
 is then used for the rest of the year, with the goal of tracking the defined
 reference temperature.
 
-Section~\ref{sec:SVGP_results} goes into details on the analysis of the Learning
+Section~\ref{sec:SVGP_results} goes into details on the analysis of the learning
 scheme using a \acrshort{svgp} Model. In this scenario, the model is first
 trained on the first five days of data, and updates every day at midnight with
 the new information gathered from closed-loop operation.
@@ -19,7 +19,7 @@ the new information gathered from closed-loop operation.
 \subsection{Conventional Gaussian Processes}\label{sec:GP_results}
 
 The first simulation, to be used as a baseline comparison with the
-\acrshort{svgp} Models developed further consists of using a `static'
+\acrshort{svgp} Models developed further, consists of using a `static'
 \acrshort{gp} model trained on five days worth of experimental data. This model
 is then employed for the rest of the year.
 
@@ -94,7 +94,7 @@ $\degree$C respectively.
 \end{figure}
 
 This large difference of performance could be explained by the change in outside
-weather (Solar Irradiance and Outside Temperature --- the exogenous inputs) from
+weather (solar irradiance and outside temperature --- the exogenous inputs) from
 the one present during the training phase. It can be seen in
 Figure~\ref{fig:Dataset_outside_temperature} that already at 500 points in the
 simulation both the GHI and the Outside Temperature are outside of the training
@@ -143,11 +143,11 @@ at midnight using the newly accumulated data from closed-loop operation.
 The results of this setup are presented in
 Figure~\ref{fig:SVGP_fullyear_simulation}. It can already be seen that this
 setup performs much better than the initial one. The only large deviations from
-the reference temperature are due to cold --- when the \acrshort{hvac}'s limited
-heat capacity is unable to maintain the proper temperature. Additionnaly, the
-\acrshort{svgp} controller takes around 250-300ms of computation time for each
-simulation time, decreasing the computational cost of the original \acrshort{gp}
-by a factor of six.
+the reference temperature are due to cold weather, when the \acrshort{hvac}'s
+limited heat capacity is unable to maintain the proper temperature.
+Additionnaly, the \acrshort{svgp} controller takes around 250 - 300ms of
+computation time for each simulation time, decreasing the computational cost of
+the original \acrshort{gp} by a factor of six.
 
 
 
@@ -161,7 +161,7 @@ by a factor of six.
 
 \clearpage
 
-Comparing the Absolute Error of the Measured vs Reference temperature for the
+Comparing the absolute error of the measured vs reference temperature for the
 duration of the experiment (cf. Figure~\ref{fig:SVGP_fullyear_abserr}) with the
 one of the original experiment, the average absolute error is reduced from 1.33
 $\degree$C to only 0.05 $\degree$C, with the majority of the values being lower
@@ -210,7 +210,7 @@ behaviour of the plant over all the experimental steps in the first two cases.
 It still has a noticeable error when predicting the behaviour of the plant on
 new data (i.e. simulations starting at steps 10750 and 11000), but it is much
 less than before. This gives a hint at the fact that the \acrshort{svgp} model's
-performance ameliorates throughout the year, but it does require much more data
+performance improves throughout the year, but it does require much more data
 than the classical \acrshort{gp} model to capture the building dynamics.
 
 \begin{figure}[ht]
@@ -361,10 +361,10 @@ Analyzing the results of a simulation done on only one day's worth of initial
 simulation data (cf. Figures~\ref{fig:SVGP_96pts_fullyear_simulation}
 and~\ref{fig:SVGP_96pts_abserr}) it is very notable that the model performs
 almost identically to the one identified in the previous sections. This
-nightlights one of the practical benefits of the \acrshort{svgp} implementations
-compared to the classical \acrlong{gp} -- it is possible to start with a more
-rough controller trained on less data and refine it over time, reducing the need
-for cumbersome and potentially costly initial experiments for gathering data.
+highlights one of the practical benefits of the \acrshort{svgp} implementations
+compared to the classical \acrlong{gp} -- it is possible to start with a rougher
+controller trained on less data and refine it over time, reducing the need for
+cumbersome and potentially costly initial experiments for gathering data.
 
 \begin{figure}[ht]
     \centering
@@ -400,7 +400,7 @@ based on closed-loop operation.
 \end{figure}
 
 As it can be seen in Figure~\ref{fig:SVGP_480window_fullyear_simulation}, this
-model is unable to properly track the reference temperature. In fact, five days
+model is unable to exhaustively track the reference temperature. In fact, five days
 after the identification, the model forgets all the initial data and becomes
 unstable. This instability then generates enough excitation of the plant for the
 model, to again learn its behaviour. This cycle repeats every five days, when the
@@ -470,7 +470,7 @@ models can be deployed with less explicit identification data, but they will
 continue to improve over the course of the year, as the building passes through
 different regions of the state space and more data is collected.
 
-These results do not, however, discredit the use of \acrlong{gp} for employment
+However, these results do not discredit the use of \acrlong{gp} for employment
 in a multi-seasonal situation. As shown before, given the same amount of data
 and ignoring the computational cost, they perform better than the alternative
 \acrshort{svgp} models. The bad initial performance could be mitigated by

90_Conclusion.tex

@@ -0,0 +1,85 @@
+\section{Conclusion}~\label{sec:conclusion}
+
+The aim of this project was to analyse the performance of \acrshort{gp} based
+controllers for use in longer lasting implementations, where differences in
+building behaviour become important compared to the initially available data.
+
+{\color{red}
+First, the performance of a classical \acrshort{gp} model trained on 5 days
+worth of experimental data was analysed. This model turned out to be unable to
+correctly extrapolate building behaviour as the weather changed throughout the
+year.
+}
+
+First, the performance of a classical \acrshort{gp} model trained on five days
+worth of experimental data was analysed. Initially, this model performed very
+well both in one step ahead prediction and multi-step ahead simulation over new,
+unseen, data. With the change in weather, however, the model shifted  
+
+Several \acrshort{svgp} implementations were then analysed. They turned out to
+provide important benefits over the classical models, such as the ability to
+easily scale when new data is being added and the much reduced computational
+effort required. They do however present some downsides, namely increasing the
+number of hyperparameters by having to choose the number of inducing locations,
+as well as performing worse than then classical \acrshort{gp} implementation
+given the same amount of data.
+
+%Finally, the possible improvements to the current implementations have been
+%addressed, noting that classical \acrshort{gp} implementations could also be
+%adapted to the \textit{learning control} paradigm, even if their implementation
+%could turn out to be much more involved and more computationally expensive than
+%the \acrshort{svgp} alternative.
+
+
+\subsection{Further Research}~\label{sec:further_research}
+
+Section~\ref{sec:results} has presented and compared the results of a full-year
+simulation for a classical \acrshort{gp} model, as well as a few incarnations of
+\acrshort{svgp} models. The results show that the \acrshort{svgp} have much
+better performance, mainly due to the possibility of updating the model
+throughout the year. The \acrshort{svgp} models also present a computational
+cost advantage both in training and in evaluation, due to several approximations
+shown in Section~\ref{sec:gaussian_processes}.
+
+Focusing on the \acrlong{gp} models, there could be several ways of improving
+its performance, as noted previously: a more varied identification dataset and
+smart update of a fixed-size data dictionary according to information gain,
+could mitigate the present problems.
+
+Using a Sparse \acrshort{gp} without replacing the maximum log likelihood
+with the \acrshort{elbo} could improve performance of the \acrshort{gp} model at
+the expense of training time.
+
+An additional change that could be made is inclusion of the most amount of prior
+information possible through setting a more refined kernel, as well as adding
+prior information on all the model hyperparameters when available. This approach
+however goes against the `spirit' of black-box approaches, since significant
+insight into the physics of the plant is required in order to properly model and
+implement this information.
+
+On the \acrshort{svgp} side, several changes could also be proposed, which were
+not properly addressed in this work. First, the size of the inducing dataset was
+chosen experimentally until it was found to accurately reproduce the manually
+collected experimental data. In order to better use the available computational
+resources, this value could be found programmatically in a way to minimize
+evaluation time, while still providing good performance. Another possibility is
+the periodic re-evaluation of this value when new data comes in, since as more
+and more data is collected the model becomes more complex, and in general more
+inducing locations could be necessary to properly reproduce the training data.
+
+Finally, none of the presented controllers take into account occupancy rates or
+adapt to possible changes in the real building, such as adding or removing
+furniture, deteriorating insulation and so on. The presented update methods only
+deals with adding information on behaviour in different state space regions, i.e
+\textit{learning}. Additionally, their ability to \textit{adapt} to changes in
+the actual plant's behaviour should be further addressed.
+
+\section*{Acknowledgements}
+
+I would like to thank Manuel Koch for the great help provided during the course
+of the project starting from the basics on CARNOT modelling, to helping me
+better compare the performance of different controllers, as well as Prof.\ Colin
+Jones, whose insights were always very guiding, while still allowing me to
+discover everything on my own.
+
+\clearpage

99A_GP_hyperparameters_validation.tex

@@ -43,14 +43,14 @@
 \begin{figure}[ht]
     \centering
     \includegraphics[width = \textwidth]{Plots/GP_313_training_performance.pdf}
-    \caption{Prediction performance of \model{1}{1}{3} on the training dataset}
+    \caption{Prediction performance of \model{3}{1}{3} on the training dataset}
     \label{fig:GP_313_train_validation}
 \end{figure}
 
 \begin{figure}[ht]
     \centering
     \includegraphics[width = \textwidth]{Plots/GP_313_test_performance.pdf}
-    \caption{Prediction performance of \model{1}{1}{3} on the test dataset}
+    \caption{Prediction performance of \model{3}{1}{3} on the test dataset}
     \label{fig:GP_313_test_validation}
 \end{figure}
 

main.tex

@@ -98,6 +98,7 @@ temperature control}
 \usepackage{fancyhdr}
 \pagestyle{fancy}
 \setlength\headheight{35pt}
+\setlength\footskip{13.6pt}
 \fancyhf{Multi-season GP performance for buildings}
 \rhead{\includegraphics[width=2cm]{Logo-EPFL.png}}
 \lhead{}
@@ -127,16 +128,14 @@ temperature control}
 \clearpage
 
 \input{10_Introduction.tex}
-\input{20_Previous_Research.tex}
 \input{30_Gaussian_Processes_Background.tex}
 \input{40_CARNOT_model.tex}
 \input{50_Choice_of_Hyperparameters.tex}
 \input{60_The_MPC_Problem.tex}
 \input{70_Implementation.tex}
 \input{80_Results.tex}
-\input{90_Further_Research.tex}
-\input{95_Conclusion.tex}
-\printbibliography
+\input{90_Conclusion.tex}
+\printbibliography[heading=bibintoc]
 \appendix
 \input{99A_GP_hyperparameters_validation.tex}
 \input{99C_hyperparameters_results.tex}