<jupyter_start><jupyter_text>Lab 6 - Gaussian Process In general, Gaussian Process (GP) is a stochastic process used for modeling data, which were observed over time, space or both. Main thing that can characterise GP is that is a kind of generalization of normal probability distributions, where each of them describes a random variable (scalar or vector if we deal with multivariate distribution).Following definition captures the essence of Gaussian Process: A Gaussian process is a collection of random variables, any Gaussian process finite number of which have a joint Gaussian distribution.If we define mean function $m(x)$ and the covariance function mean function $k(x,x')$ of a real process $f(x)$ as:\begin{equation}m(x) = E[f(x)]\end{equation}\begin{equation}k(x,x') = E[f(x)-m(x))(f(x')-m(x'))]\end{equation}Then we can define GP can as:\begin{equation}\label{eq:gpmain}f(x) \sim GP(m(x),k(x,x'))\end{equation}GP can be fully determined by only declaring mean and covariance functions. Mean in most cases is set to value "0", because such a setting can be useful, simplifies matters and is not a difficult requirement to fulfill. Of course, there are exampleswhere we would like to change the mean e.g. for better model interpretability or the specification of our prior. Covariance function (also called kernel function) representsa similarity between data points. Note that usually covariance is chosen fromthe set of already defined functions. The one should at least pick one, which representsprior beliefs of the problem. But in fact, covariance function can be any function, whichhas the property of generating a positive definite covariance matrix. Anyway, creating and defining new covariance functions, which will simultaneously be correct and have apractical usage, can be really difficultThe most basic and common kernel function is Radial Basis Function (RBF), which is defined with the formula:\begin{equation}\label{eq:rbf}k\left(x_{i}, x_{j}\right)=\exp \left(-\frac{d\left(x_{i}, x_{j}\right)^{2}}{2 l^{2}}\right)\end{equation}Its main property is that its value is usually only dependent on the distance from the specified point. The parameter which RBF kernel uses is $l$ as a characteristic length scale. RBF is infinitely differentiable. That means the GP with this kernel has a mean square derivatives for all orders. Other worth noting kernel, in terms of our classess is one from exponentiated quadratic family, where $\rho$ is length-scale and $\alpha$ is marginal deviation:$$k\left(x_{1}, x_{2}\right)=\alpha^{2} \exp \left(-\frac{1}{2}\left(\frac{\left|x_{1}-x_{2}\right|}{\rho}\right)^{2}\right)$$Gaussian Process is growing in popularity of usage in domains connected with statistics and machine learning. It is especially efficient solution, while we are dealing with a problem of obtaining data e.g. for faulty cases or just an overall lack of comprehensive data and if we want to obtain not only prediction but also confidence interval. Main areas where GP is used for are: regression, prediction, classification and identification. On this classes we will focuse on a simple task of regression.  Task 1: simulating from Gaussian Process model To make any simulations, first you need a set of observations:<jupyter_code>#Observed data
x_obs = [-10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10]
y_obs = [0.328572824089476, 4.20607004111644, 1.35507551134795,
0.161608755204364, -5.42320349780782, -3.05851276224202, -0.0764172642034502,
-4.55218472276499, -0.902226297922731, -5.8609833528976, -1.05854090910473]<jupyter_output><empty_output><jupyter_text>Model also needs sampling space and the indexes of where our observed data resides there. For your convenience, we provide them below:<jupyter_code>x_predict = np.linspace(-11,11,551,endpoint=True)
observed_idx = [26, 76, 126, 176, 226, 276, 326, 376, 426, 476, 526]<jupyter_output><empty_output><jupyter_text>In this example we are using kernel from exponentiated quadratic family, so we need to provide its hyperparameters (alpha and rho). Moreover, for sample generation, we would like to set a sigma value. Let's start with some arbitrarily selected parameters:<jupyter_code>alpha = 3
rho = 5.5
sigma = 2<jupyter_output><empty_output><jupyter_text>Now we can create a model and sample from it. Do not forget to define proper data value based on requirements from stan model.<jupyter_code>model = CmdStanModel(stan_file='gaussian_process.stan')
result = model.sample(data=data, chains = 4)<jupyter_output><empty_output><jupyter_text>Your task now is to:- plot few (about 5-10) samples generated from model- plot mean value from model with its confidence interval (standar deviation, using errorbar function)Also on each plot place points with obsserved data. Repeat the steps for two other parameters sets of your choice  Task 2: optimize hyperparameters To achive better results of GP fitting for our observed data, we need to oprimize its parameters. We can do that by using maximum marginal likelihood estimation. It can be done by placing parameters in "parameter" block in stan model and using "optimize" method.<jupyter_code>model2 = CmdStanModel(stan_file='gaussian_process_optimize.stan')
result2 = model2.optimize(data=data, seed=5838298)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.08045540796963947, 'fraction_numerical': 0.04800759013282732, 'mean_word_length': 5.630188679245283, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 21, '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': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '610042', 'license': 'CC-BY-4.0', 'n_tokens_mistral': 1605, 'n_tokens_neox': 1402, 'n_words': 673, 'path': 'Data Analytics 2022/Lab 6 - Gaussian Process/Lab6.ipynb', 'repo_name': 'piotrbullet/public_lectures', 'size': '5270'}